VERTICAL CURRENT STRUCTURE IN THE GREAT LAKES
Vincent E. Noble
Joseph C. Huang
James H. Say lor
Final Report
Federal Water Pollution Control Administration
Grant WPOIO67
Special Report No, 57
of the
Great Lakes Research Division
The University of Michigan
Ann Arbor^ Michigan
1968
TABLE OF CONTENTS
Page
LIST OF TABLES v
LIST OF FIGURES vi
INTRODUCTION 1
I. CHARACTERIZATION OF THE CIRCULATION DYNAMICS OF LAKE MICHIGAN
(Vincent E. Noble) k
II. CURRENT METER STUDIES OFF SHEBOYGAN, WISCONSIN, AUGUST 196?
(Vincent E. Noble) 12
III. DETERMINATION OF EDDY VISCOSITY AND EDDY DIFFUSIVITY IN LAKE
MICHIGAN (Joseph C. Huang) 26
Introduction 26
The Magnitude of Eddy Viscosity and Eddy Diffusivity in Lake
Michigan 28
An estimate of eddy viscosity by Reynolds stress and vertical
velocity gradients 28
An estimate of the upper bound value of eddy viscosity by the
result from statistic theory 29
An estimation of vertical eddy viscosity from current
measurement near the bottom 29
An approximation of the vertical eddy viscosity by surface
wind 51
A measurement of the eddy diffusivity by the dispersion of
marked particles 3^
Instantaneous sources 55
Continuous sources 55
A measurement of the horizontal diffusion by drogue studies 58
Comparison of the Values of Eddy Viscosity and Eddy Diffusiv
ity with Other Studies 58
Conclusion ij5
IV. NUMERICAL EVALUATION OF STERN'S CIRCULATION MODEL (James. H. Say
lor) k6
Development of Stern *s Model k6
Application of Stern's Model to Lake Michigan 5I
Effect of density difference, Ap 55
Effect of T 56
Effect of depth of thermocline 56
iii
TABLE OF CONTENTS (Concluded)
Page
V. THERMAL CURRENT STRUCTURE IN LAKE MICHIGAI^, ASSOCIATED WITH THE
SPRING WARMUP SEASON— A THEORETICAL STUDY (Joseph C. Huang) 59
Introduction 59
Modeling Considerations and Assumptions 6o
Mathematical Formulation 6k
Nondimensional equations 65
Solutions of the Zeroth p Order Approximation 68
Temperature 68
Boundary layer scaling 68
Interior solutions 69
Boundarylayer solutions 71
Results of the Zeroth pOrder Solution 7^
Temperature 7^
Velocity 85
The First p Order Solutions 86
Conclusions 91
Future Plan 92
IV
LIST OF TABLES
Tables Page
51 « Eddy viscosity on two currentmeter buoy stations (172^^ August
1967 at 8 miles off Sheboygan). 50
52 « Estimated vertical eddy viscosity by current measurement near
bottom, 51
55* Estimated vertical viscosity by wind speed, 52
54„ General information and effective eddy diffusivity of dye patch
studies. 5U
55« General information and effective eddy diffusivity of continu
ous source dye plume study. 59
56. General information and eddy diffusivity of drogues, study in
Lake Michigan. k2
h"!^ Calculated eddy diameters, phase velocities, and internal wave
periods of upper layer in Lake Michigan. 55
U2. Calculated eddy diameters, phase velocities, and internal wave
periods for upper layer in Lake Michigan with density difference
changed. 55
^5< Calculated eddy diameters, phase velocities, and internal wave
periods for upper layer with wind stress changed. 56
kk. Calculated eddy diameters, phase velocities, and internal wave
periods in Lake Michigan with the depth of thermocline changed. 56
51. Nomenclatures and nondimensional parameters. 65
52. Local coordinates in different regions. 70
V
LIST OF FIGURES
Figure Page
11. Surface temperature contours from Airborne Radiation Thermometer
flight of 18 October I966. 5
12. Surface temperatures and 010 m dynamic height currents from BT
transects taken by R/v INLAND SEAS, October I966. 6
15. Surface temperature contours from Airborne Radiation Thermometer
flight of 25 October 1967. 8
1U. Surface temperatures and 010 m dynamic height currents from BT
transects taken by R/V INLAND SEAS^ 26 October I967. 9
15, Temperature structure of a vertical crosssection of Lake Michi
gan corresponding to Airborne Radiation Thermometer flight
tracks. 10
2la Location of current meter stations, showing contours of bottom
topography. I5
22. Comparison of Milwaukee (MKE) winds and il/2 hour current vec
tor averages for stations I and II at 010 m depths (15 August 
2k August 1967). 15
25. Dynamic height current components and corresponding ^l/j hour
current meter averages, I7 August I968. I7
2U. Drogue observations and quarterperiod current meter averages,
18 August 1968. 18
25. Drogue observations (showing drogue displacement and net trans
port velocity), dynamic height current components, and corre
sponding ^1/5 hour current meter averages. I9
26. Dynamic height current components and corresponding ^l/5 hour
current meter averages^ 21 August I967. 21
27. Dynamic height current components and ^l/5 hour current meter
averages, 23 August I967. 22
51. A typical distribution of dye concentration across a dye patch. ^k
VI
LIST OF FIGURES (Concluded)
Figure Page
52. Lateral growth of the dye patch with time. 55
55* Longitudinal growth of the dye patch with time. 5^
5^. Continuous source apparatus. 57
55* Lateral growth of dye plume at a constant distance from the
source. 40
56. Lateral growth of dye plume with distance from source. kl
i+l. Schematic twolayer model. 51
U2. Specific parameters of the twolayer model. 5^
51. Eastwest cross section of the central portion of Lake Michigan. 6l
52. Cross section of a long, symmetrical trapezoidal Lake Michigan
model. 62
55 Boundary regions of the Lake Michigan model. 70
5^. Surface temperature patterns during spring warming of Lake
Michigan from available data, 195^68. 75
55^ Vertical temperature structure of central portion of Lake
Michigan during warming period, I96568. 79
56. Generalized horizontal zonal currents. 85
57» Meridional and vertical circulation. 86
5^8. Surface currents in Lake Michigan, 29 June 1955 (after Ayers
et al. 1958). 87
59» Brogue studies for surface current measurements near eastern
shore . 88
510. Drogue studies for surface current measurements near western
shore . 90
511. Circulation associated with "Thermal Bar" (after Rodgers I965). 91
VI 1
INTRODUCTION
The research program carried out under this research grant (WPOIO67) was
a twoyear continuation of the studies initiated under a previous U.S. Public
Health Service grant (WP0079^). The final report describing the previous work
was prepared by Noble (I966). The previous report shows the failure of Ekman
type models in synthesizing current meter records from buoy stations when using
the wind measured at the station as an input stress. The previous study further
suggests that geostrophic forces defined by the temperaturedependent density
structure of Lake Michigan may be the dominant factor in determining the three
dimensional circulation patterns of the lake.
Four major field experiments were carried out (in cooperation with other
projects and individuals) under this grant. The two most significant operations
were thermal mapping flights carried out by the Ant i Submarine Warfare Environ
mental Prediction Services (ASWEPS) of the U.S. Naval Oceanographic Office. The
results of the l821 October I966 flight program of the ASWEPS Super Constella
tion^ EL COYOTE (showing the existence of a major circulation feature^ the
Coyote Current) have been described in detail in a report by Nobel^ £t al. (1967).
The results of the second flighty on 25 October 1967^ are compared with those
of the 1966 flight in a paper in preparation by V. E. Noble and J. C. Wilkerson.
A third major field experiment was carried out in August 1967j> in coopera
tion with Dr. G. E. Birchfield of Northwestern University. A summary of the
results of this study was presented by Birchfield and Noble. ^ In this study^
two current meter stations were established in Lake Michigan for a twoweek
period^ and comparisons were made between the current meter records, drogue
measurements, and dynamic height current calculations. Spectral analyses and
inertial period averages were also computed from the current meter records.
The fourth field experiment consisted of a series of measurements of the
temperature, eddy viscosity, eddy diffusivity, and current structures in the
nearshore edge of the lake during the spring warming period of I968. These
measurements were carried out by J. C. Huang. The results will be used as part
of the basic input data and parameters for his theoretical investigations of
lake circulation, and will be compared with the results obtained in connection
with another project by Noble and Anderson (1968) during the spring of I967.
In addition to the experimental program, two theoretical studies have been
carried out in connection with this grant by students of the Department of
Statistical analysis of currents at two nearby stations in Lake Michigan, summer
1967.. G. E. Birchfield and V. E. Noble. Paper presented at 11th Conf . on
Great Lakes Res., Milwaukee, I968.
Meteorology and Oceanography of The University of Michigan. The first was the
output from a Special Problems' course in which James H. Saylor calculated
numerical results, based on Lake Michigan data, from Stern's (1965) model for
the interaction of a uniform wind stress with a geostrophic vortex. The second
is a Ph.D. dissertation being written by Joseph C. Huang. Mr. Huang is doing
a theoretical analysis of a simplified Lake Michigan model using the early
spring heat input to the lake as the prime current gene rating force. Early re
sults from this model (Huang's dissertation is not expected to be completed
until the summer of I969) indicate that the geostrophic currents resulting from
the spring temperature structure of the lake may be as great as previous compu
tations for winddriven currents. Those results correspond with observations,
and with the concepts within Stern's (I965) model.
In the proposal for this grant, it was anticipated that the data from a
special triangular array of current meter stations set in Lake Huron during the
winter of I96566 would provide a basis for investigation of the Eulerian
Lagrangian transformation inherent in the extrapolation of current meter data,
and of correlations between data from pairs of current meters. Instrumentation
problems caused a lack of data from this source, and therefore it was necessary
to rely upon measurements made under this grant for subsequent analysis and in
terpretation.
References
NOBLE, VINCENT E. I966. VERTICAL CURRENT STRUCTURE IN TEE GREAT LAKES • Spec.
Rep. No. 27, Great Lakes Res. Div., University of Michigan, k2 p.
NOBLE, VINCENT E., et al. I967. NAVOCEANO FLIGHT OVER LAKE MICHIGAN, 1721
October I966. Spec. Rep. No. 5I, Great Lakes Res. Div., University of Michigan,
17 p. + 2k figs.
NOBLE, VINCENT E. and R. F. ANDERSON. I968. Temperature and current .structure
in the Grand Haven, Michigan, vicinity during thermal bar conditions. Proc.
11th Conf. on Great Lakes Res., Internat. Assoc, for Great Lakes Res. In press.
STERN, MELVIN E. I965. Interaction of a uniform wind stress with a geostrophic
vortex. DeepSea Research, 12:555567.
I. CHAEACTERIZATION OF THE CIRCULATION DYNAMICS OF LAKE MICHIGAN
Vincent E. Noble
The results of previous investigations carried out in connection with this
program tend to indicate that the circulation dynamics of Lake Michigan and, by
inference, the other Great Lakes, are dominated by geostrophic forces rather
than the wind stresses during the summer season. Analysis of temperature records
from the Great LakesIllinois River Basins Project thermographs installed at
buoy stations in Lake Michigan during the winter season of 19656^ indicated
that the prevailing westerly wind drift caused an eastward surface current on
Lake Michigan, a sinking on the east side of the lake basin, a westward bottom
current, and an upwelling on the west side of the basin (Heap and Noble I966) .
In the spring, the "thermal bar" (Rodgers I965; Noble and Anderson I968), which
is a consequence of the spring heating of the lake, defines a general circula
tion around the edge of the lake basin. During the onset of the summer season,
with the development of the thermocline, the wind stressgeostrophic stress
interactions become extremely complex.
Critical testing of theoretical models depends heavily upon our ability to
make adequate density measurements in both time and space. Analysis of
data from single current meter stations has failed to demonstrate a direct rela
tionship between wind and current. Not enough is yet known about the spatial
structure of currents to interpret the relationship between data from current
meter stations separated by distances of 5 to 10 miles. Surface temperature
patterns obtained by synoptic observations from aircraft imply certain types of
circulation features when interpreted with respect to dynamic height current
computations made from BT transects. However, this interpretation is subject
to errors inasmuch as the dynamic height computations depend on the essential
assumption that the current field is defined by the density (temperature) field,
and therefore it is only selfconsistent that the dynamic height currents should
show a strong correlation with the surface temperature structure. To the extent
that dynamic height computations show a good qualitative correlation with cur
rents measured with drogues and current meters, the fundamental assumption of
the dominance of the geostrophic forces is valid. The question yet to be re
solved is the interaction between the wind and geostrophic forces. The persis
tence of small details in the surface temperature structure of the lake for
periods of three weeks (Heap and Noble I966; Noble I967) give credence to the
postulate of dominance of the lake circulation dynamics by geostrophic forces.
Figure 11 shows isotherms contoured from 5mile averages of surface
temperature along 25 eastwest flight lines spaced at 10mile intervals from
the south end of Lake Michigan on 18 October I966. Subsequent flight operations
showed that the general temperature pattern persisted for at least three subse
quent days. Figure 12 shows the dynamic height currents computed for the 010 m
42*00'
8e»oo'
Figure 11. Sairface temperature contours from Airborne Radiation
Thermometer flight of 18 October I966.
20
10
10
I
87» 35'
&?• 15'
I
86« 55'
I
86» 35*
INLAND SEAS 18 OCTOBER 1966
43* 57.2'
87* 35'
87»I5
86» 55'
INLAND SEAS 20 OCTOBER 1966
43* 57.2*
— 44* 00*
87* 35'
87* 15'
86» 35*
INLAND SEAS 19 OCTOBER 1966
Figure 12. Surface temperatures and 010 m dynamic height currents
from BT transects taken by R/v INLAND SEAS, October I966.
layer, compared with surface temperature readings taken from BT transects along
the 45 '^ 57, 2* N flight line on 18, 19, and 20 October I966 by R/v INIAp SEAS.
These dynamic height currents indicated a strong northward current parallel to
the strong temperature gradient on the west side of the northern half of the
lake. This circulation feature was also indicated by a distinct visual water
mass separation along the gradient. The water was warm and green to the east
of the gradient, and cold and blue to the west. The circulation feature ^as
named the Coyote Current in honor of the NAVOCEANO research aircraft, the Super
Constellation EL COYOTE.
If the circulation of Lake Michigan is dominated by geostrophic (temperature
dependent density) forces, it would be expeated that the Coyote Current would
appear each fall in the northern end of Lake Michigan. Therefore another flight
was scheduled during the fall of I967 to see if the current did recur.
The early fall period of I967 was much colder than that of 1966, and the
temperature mapping flight could not be scheduled until 25 October. Therefore
the 1967 flight was equivalent to being two weeks later in the season than that
of 1966. The surface temperature contours from the I967 Coyote flight are shown
in Figure I5. Due to heavy weather, the INLAKD SEAS could not make a BT transeci
on the day of the flight. Figure lk shows 010 m dynamic height currents and
the surface temperature profile from the INLAM) SEAS BT transect at ky^^J.^^'N
on 26 October I967.
The 1967 temperature gradient was stronger (of the order of ^''C). The
temperature gradient turned offshore near Sheboygan, Wis., in I966, and more
nearly paralleled the shoreline in I967. Figure 15 shows the temperature struc
ture of a vertical section of the lake as determined by the BT transects of
1966 and 1967 Because the development of the surface temperature structure
was not followed continuously through the fall and early winter periods of 1966
and 1967^ it is not possible to make definitive judgments based on the compari
son of measurements made on 18 October I966 and 25 October I967. However, the
appearance of the strong temperature gradient feature covering a large portion
of the north south extent of the lake basin, accompanied by a water mass separa
tion as evidenced by color changes visible from the aircraft, gives a strong
indication that a similar phenomenon was observed in both years.
Additional smallscale experiments were planned and executed to provide
experimental foundations for testing and evaluation of theoretical models that
would lead to an understanding of the circulation dynamics of the lake. These
experiments were continued observations of the current regime associated with
the early spring temperature structure, and a 2week period of measurement using
two current meter stations and comparison current drogue and BT observations.
Figure 15 • Surface temperature contoiirs from Airborne Radiation
Thermometer flight of 25 October I967.
10
1 11 'I '■ '
10
43«'57.3'
INLAND SEAS
26 OCTOBER 1967
Figure 1k. Surface temperatures and 010 m dynamic height currents
from BT transects taken by r/v INLAND SEAS, 26 October 196?.
(a) INLAND SEAS 18 OCTOBER 1966
go 70 QO 90 IJO
50 5^<C^
100
150
==;^
\ ^^
7=^ ';'C^^ 2 ^'^'^^^u***^''
87«30'
I
8 7" 00'
I
(b) INLAND SEAS 20 OCTOBER 1966
go go 90 0
86<'30
I
43''57.3
(c) INLAND SEAS 26 OCTOBER 1967
503070 e^S* , 10" 9" 9°
87*'30'
87*00'
86<»30
Figure 15. Temperature structure of a vertical crosssection of Lake
Michigan corresponding to Airborne Radiation Thermometer flight tracks,
(a) 18 October 1966. (b) 20 October 1966. (c) 26 October I967.
10
References
HEAP, J. A. and V E. NOBLE. I966. GROWTH OF ICE ON LAKE MICHIGAN. Spec.
Rep. No. 26, Great Lakes Res. Div., University of Michigan, 9^ p.
NOBLE, V. E. 1967 Evidences of geostrophically defined circulation in Lake
Michigan. Proc. 10th Conf. Great Lakes Res., Internat. Assoc, for Great Lakes
Res., p. 289298.
NOBLE, V. E. and R. F. ANDERSON. I968. Temperature and current structure in
the Grand Haven, Michigan, vicinity during thermal bar conditions. Proc. 11th
Conf. on Great Lakes Res., Internat. Assoc, for Great Lakes Res. Jn press.
RODGERS, G. K. I965. The thermal bar in the Laurentian Great Lakes. Univer
sity of Michigan, Great Lakes Res. Div. Pub. No. 15:558563.
11
II. CURRENT METER STUDIES OFF SHEBOYGAN, WISCONSIN, AUGUST I967
Vincent E, Noble
During the period of I6 August through 2k August IS&J , two current meter
buoy stations were established in 280 ft of water about 8 miles offshore from
Sheboygan, Wis. The station locations and relative bottom topography in that
area of the lake are shown in Figure 21. The buoy stations were set a half
mile apart in roughly an eastwest direction with current meters set nominally
at the surface and at 10m depths. The actual current speed rotors were at 2
and 12 m, respectively. Both current meter stations were obtained on loan from
the Federal Water Pollution Control Administration, Great Lakes Illinois River
Basins Project. The experiment was designed to compare current drogue obser
vations and dynamic height current computations with current meter readings, and
to provide an insight into the distance through which current meter readings
can be extrapolated and into the EulerianLagrangian transformation represented
by current drogue and current meter measurements.
Station I was set at l^ii^O EST on 16 August I967 at a position of kykh.VE
and 87''32.0'W. The station depth was 2&\ ft. Station II was set at I530 on
16 August 1967 at a position hyhk.2''^ and 87^52. 6'W. The water depth at sta
tion II was 280 ft. The station position was onehalf mile at bearing 287*^ from
station I. The movie film on which the data from the surface meter at station
I were recorded was prematurely exposed, therefore, no useful data were obtained
for the first ^15 hours after the station was installed. The other three current
meters operated satisfactorily throughout the full period.
The buoy station configuration was similar to that used by the Great Lakes
Illinois River Basins Project. In this instance, however, the subsurface buoy
(which is used to provide tension on the meter string) was set so that it was
just awash at the surface of the water. The thermocline was rather shallow in
the region of measurement so that the 12m current rotor was just about at the
thermocline, generally being in the metalimnion. The weather during the period
of measurement was generally unfavorable so that far less comparison data were
obtained than had been desired. On I7 August only BT's were taken for dynamic
height current computations. On the l8th only drogue observations were made.
On the 19th the ship was weatheredin at the dock. On the 20th both drogue and
BT observations were made. On the 21st, BT*s only. The ship was laid in on the
22d, with only BT*s again on the 25d. The stations were recovered on the 2^th.
The winds were generally southwest to west on the l6th and 17th, with a frontal
passage on the l8th producing northwest winds through the 19th. The wind went
back to southwest to west on the 20th with a second cold front passage with
north to northeast winds on the 21st remaining northeast through the 22d and
going southeast on the 25d.
12
4^50
87'50'
SHEBOYGAN
■43M0
Figure 21. Location of current meter stations, showing
contours of bottom topography.
15
The full statistical analyses of the current meter data have been discussed
by Birchfield and Noble, The statistical analyses^ performed by Dr. G. E*
Birchfield, Department of Engineering Sciences, Northwestern University, include
determination of the mean transport current and a measure of the dispersion of
the instantaneous current vectors during full and quarter inert ial periods (17I/5
hours, and ^l/5 hours, respectively). In this report, the quarterperiod means
are compared with drogue measurements, dynamic height current computations, and
wind records from Milwaukee and Muskegon. The wind records from the land stations
are used to document the general wind regime and significant frontal passages.
Figure 22 compares the wind records from the land stations with the quarter
period averages from each of the current meters.
The quarterperiod current meter averages discussed in this section of the
report are the result of a first analysis carried out by Dr. Birchfield. The
quarterperiod averages were computed by taking the times that the current meters
were set as t = 0. Therefore, the 4l/5hour averages for stations I and II
are displaced in time by 50 minutes (the elapsed time between the setting of
stations I and II) . The current meter averages are being calculated for coin
cident periods, and a complete report of the statistical analysis of the current
meter data and drogue comparisons is being prepared by Birchfield and Noble.
The complete report will discuss the current meter averages and vector disper
sions for the coincident periods, and will contain spectral analyses of the
current meter data.
The winds at Milwaukee (MKE) were S from I5 August through 0000 I8 August.
A frontal passage occurred around O5OO on I8 August, with winds shifting to N
at 15 kts (ships on the lake reported 1825 kts, N) . The MKE winds diminished
on the 19th and went light, SW on the 20th (a ship reported JO i^"ts at 190° at
2100 on 20 August). MKE winds were WSW, 510 kts through 1200 21 August, with
a strong frontal passage around I5OO on 21 August (ships reported 1555 kts at
060°020° behind the front during the period 09002100 on 21 August). MKE
winds went east, decreasing through I80O on 2k August. (Ship reports gave 21 kts
at 030° for 2100 on 22 August, and 09 kts at 010° for 2100 on 25 August.) The
ship reports used for this study were carferry reports taken approximately 50
miles NE of the buoy stations.
As shown in Figure 22, the quasi steady currents between the fronts do not
show a simple relation to the wind direction. The response of the current meter
records to a frontal passage is the development of a rotational characteristic
lasting for a period of I8 hours. (The rotational characteristic develops after
a time lag of the order of 12 hours.) It is suggested that the complex relation .
between the current vectors and wind field may be explained in Ekman transport
and acceleration of geostrophic eddies as described by Stern ^s (I965) model.
2
See footnote 1, p. 1.
Ik
Figure 22 Comparison of Milwaukee (MKE) winds and 41/2 hour current vector
averages for stations I and II at 010 m depths (15 August  24 August 1967)
Figures 25 through 27 show the results of the comparative measurements
between the drogues, dynamic height current computations, and the quarter iner
tial period averages of the cu!i:rent meter records. The quarter period averages
give the average current vector determined from a 4l/5hour interval. The
dynamic height current computations yield currents between pairs of BT*s taken
within a time interval of the order of l/2 hour within the comparison current
meter period. The net drogue movements indicated by the dotted lines on the
figures are used to determine a net transport vector.
Figure 25 shows the relation between dynamic height currents and quarter
period averages of the current meter records for 17 August 1967* For station
I the averaging interval was O80O to 1220 hours and for station II, 085O to
1510 on 17 August (period 5) The surface current data for station I were lost
due to premature exposure of the film but records are available for 10 m at
station I, and and 10 m at station II. The average currents from stations I
and II are 2O5I cm/sec. Dynamic height currents are computed from three BT
casts taken between 0850 and 0952 on I7 August. The dynamic height computations
give only the component of current speeds normal to the path between BT casts.
For this comparison period the components determined by dynamic heights were
k.9 and 5*6 cm/sec, respectively, for the average current in the to 10m
level between the three BT casts.
The dynamic height currents were calculated by the method of Ayers (1957),
and assume a reference depth of 70 m. The computed speed components were within
20^ of those shown by the average values for the current meter stations. The
dynamic height component computed for a position 2 miles east of the stations
showed good agreement with the current meter average, but the component computed
for a position II/2 miles west of the stations showed a direction I80** from
that given by the station data. The ship winds at 0815 on I7 August were I9
kts, with direction 195 °«
Figure 2k shows a comparison between drogue observations and quarter period
current meter averages for 18 August 1967. The four current drogues were set
between 0955 and IO55 and recovered between 1^^15 and 15^11 on I8 August. A frontal
passage had occurred Just before setting the drogues, and only two positions of
the drogues were determined — the set point and the recovery point. Drogue speeds
were generally to the north with a net transport velocity ranging from 57 to
21.5 cm/sec. For station I a quarter period average encompassed the time be
tween 1000 and 1^20 on I8 August with surface current of 19*5 cm/sec and 10m
current of 55*6 cm/sec (period 11). The surface current speed was 5O.8 cm/sec
with 10m speed of 52.8 cm/sec at station II.
The current drogues consisted of ^ft square sheet metal current crosses
suspended from modified net stake floats. The center of the current cross was
at about the same depth as the current speed rotor for the surface meters. The
drogue showed a good correlation with both speed and direction with the surface
current meter of station II but very poor correlation with the direction of
station I. The winds were generally north and at lUoO the ship wind was north
16
^
^
GO
00
a>
0)
c>
^
o
0)
E
K
3
O
Z
1
1
UJ
UJ
(o
cc:
o:
o
H
LlI
UJ
<
X
K
h
J
o
UJ
UJ
CL
S
:e
CO
UJ
X
o
o
z
u
o
LiJ
J
Z
z
UJ
UJ
o
UJ
UJ
3
3
:s
q:
q:
O
o
<
cr
q:
o
o
z ,
3
3
q:
cr
>" H
O
O
o
o
O GO
^ u
o
CVJ
to
00
bD
•\
a
w
•H
■p
rd
^
a
o
c^
ft
H
CQ
TiJ
^1
fl
U
•H
o
:S
o
'd
^
a
00
ai
On
CQ
H
P
a
P
(D
GQ
o
gj
f
<
o
o
H
•^
a
^
(D
CQ
^
cu
^
bD
,
o
^
o
P
a)
p
\>
CQ
^
5
"^^v..^
bD
S
•H
Sh
o
(D
^
P
a
^
•H
O
^
4^
d
a
SU
(U
0)
ft
S
^
^
CQ
p
o
a
•
0)
CVJ
g
^
o
(U
CVl
•\
S
H
o
ir\
•H
1
CA
F^
=t
H
rO
00
IT
LEGEND
O
METER  Surfoce
METER  10 meter
CURRENT
CURRENT
DROGUE
DROGUE DISPLACEMENT
DYNAMIC HEIGHTS
BT
32.8
43* 46'
19.3
43044,
1/2
87» 34'
miles
43042.
87" 30*
Figure 2k. Drogue observations and quarter period current meter
averages, 18 August I968. The net drogue displacements and average
transport velocities are shovra for a 5hour interval. Wind N, 19
25 kts.
18
^ o
I I
(T a:
ill UJ
I H
UJ UJ
UJ
o
<
J
a.
CO
<0
o
ro
o
fO
a
z
Ui
o
111
J
Z Z UJ UJ
UJ UJ 3 ID
a: q: o o
(T cr o o
3 3 q: (T
o o Q Q
2
<
z
O CQ
CM
ro
^ r t
00
CM
• e
GO
O
fO
00
O
H
I 00
CM
fO
o
00
(M
H
•H
is:
O
o
•H
(D
P
+5
CQ
o
CM
O
H
I
O
ON
O
19
CD
o
o
ro
CVJ
o
ro
0)
0)
o
i^
o
o
k
£
H
3
O
2
(f)
\
1
lU
:^
UJ
C/)
QC
oc
o
h
UJ
UJ
<
X
h
H
J
o
UJ
UJ
a.
UJ
X
Q
Q
LlI
LU
_J
2
H
Z
LU
LU
o
liJ
LU
3
3
s
q:
q:
O
CD
<
cr
q:
o
O
z
3
3
q:
q:
>■ K
O
O
Q
Q
Q CD
ilro
CvJ
CO
c\J
■ o
P
CM
CO
TIJ
a
•H
^
6
•
ro
(U
1^
CD
o
•H
•
^d
'5'
(U
j^
0)
nd
PM
CQ
H
O
4J
CJ
d
O
<u
O
^
CVJ
^
ITv
ro
o
1
o
CM
CX)
ON
(D
H
P«4
P
CQ
O
OJ
o
^
H
ro
On
o
H
00
S
I r t
20
LEGEND
CURRENT METER  Surface
CURRENT METER  10 meter
T5 DROGUE
, DROGUE DISPLACEMENT
> DYNAMIC HEIGHTS
O BT
o
43''46'
2.7
18.3
3.9
0\.^
2.4 _ o
87°34'
1/2
— miles
O
87° 32'
43° 42'
87'>30'
Figure 26. Dynamic height current components and corresponding
i^l/5 hour current meter averages, 21 August I967. Current speed
in cm/sec. Wind S, 7 kts.
21
LEGEND
CURRENT METER  Surface
CURRENT METER  10 meter
T' DROGUE
— X DROGUE DISPLACEMENT
a* DYNAMIC HEIGHTS
O BT
51.4
/
25.9 21.3
4.2
\=
43°46'
O
43044.
33.7
1/2
• miles
86° 34'
43° 42'
86°30'
Figure 27. Dynamic height cairrent components and UI/5 hour current
meter averages, 23 August I967. The dynamic height coraponents were
computed for a time at the end of averaging period 57 (solid vectors)
and the beginning of period 58 (broken vectors). Wind SE, 15 kts.
22
1925 kts. It is to be noted that all four drogues, and the surface current
averaged at station II for the quarter period, were all to the north, that is,
into the wind.
Figures 25^ and b show a comparison between drogues, dynamic height cur
rents, and current meter averages for 20 August 1967* The five drogues were set
on an eastwest line through the current meter stations between 0912 and 0935*
The drogues were fixed between 1525 and l4^0 and recovered between I655 sind 17^0.
Figure 2 5a shows the first half of the measurements between the set and
the first fix of the drogues, and Figure 25b shows the second half between the
fix and recovery. Figure 2 5a gives the quarter period average currents for
the time period between 09^0 and 1^00 for station I (period 22) and between
1050 and 1^50 for station II (period 22). The dynamic height current component
for a position 1 mile west of the current meter station shows good agreement
(speed: 0.7 cm/sec) for the northerly component with the net transport of the
drogue set at the same location (8.6 cm/sec to the north of east). The speed
component (to the north at 0.1 cm/sec) computed for a position II/2 miles to
the east of the stations is in the proper sense, but does not agree as well
with its corresponding net drogue transport (9.8 cm/sec to the northeast). Of
particular interest is the net south transport of the drogue set 2 miles to the
west of the stations. There is good agreement between the drogues set at the
current meter sites and the ^l/jhour meter averages.
Figure 2 5b compares the quarter period current meter records from 1^00 to
1820 for station I (period 25) and 1^50 to I9IO for station II (period 25) with
the net transport currents given by the motion of the drogue between their fix
and recovery points. These two direct current measurements are compared with
the to 10m current calculated by the dynamic height method. In this case,
there is less correlation between the surface current meter averages, and be
tween these averages and the drogue transports and the dynamic height current
components. (Only two dynamic height components are calculated for this inter
val, shown as 2.5 and 0.8 cm/sec to the east.)
Figure 26 gives the results of the experiments on 21 August I967 comparing
the 10m dynamic height currents with the quarter period average from station I
in the interval of 0720 to 11^0 (period 27) and the interval of 08IO to I25O
at station II (period 27). Figure 27 gives the results for 25 August in which
dynamic height currents for the to 10m layer are compared with current meter
averages for two adjacent Uhour periods. The BT casts for the computation of
the dynamic height currents were taken at the end of period 57 in the averaging
sequence, therefore, the two adjacent periods 57 and 58 are compared with the
dynamic height currents. For station I the two intervals are 02J^0 to 07OO and
0700 to 1120, and for station II O55O to 0750 and O75O to 1210 on 25 August.
Of interest in this figure are the strong differences between the dynamic height
current components and between the two sets of ^hour averages that bracket the
time interval when the dynamic height components were computed.
25
The results of these experiments, notably the results of l8 and 20 August,
raise strong questions about the extrapolation of current data from any given
single current meter station. There are often serious discrepancies between
the ^hour averages in corresponding periods between the two current meter sta
tions that are only onehalf mile apart. Further, the drogues move near the
current meter stations and show net transports over comparable ^hour periods
that are at variance with the average currents computed at the stations.
Finally, the sequence of quarterperiod current meter averages shown for
the cases of the frontal passages on l8 and 21 August show that the dynamics of
the response of currents to winds are not well understood. Moreover, in the
period preceding the front of 21 August when the winds were from the southerly
direction, the general tendency of the current meter data was to show an overall
southward drift.
2k
References
AYERS, J. C. and ROGER BACHMNN. 1957 Simplified computations for the dynamic
height method of current determination in lakes. Limnol. Oceanog. 2:155157,
STERN, MELVIN E. 1965* Interaction of a uniform wind stress with a geostrophic
vortex. DeepSea Research 12:555367.
25
III. DETERMINATION OF EDDY VISCOSITY AND EDDY DIFFUSIVITY IN LAKE MICHIGAN
Joseph C. Huang
INTRODUCTION
Theoretically, a flow field of incompressible Newtonian fluid can be
solved by the equation of conservation of mass and by the NavierStokes equa
tion of motion, satisfying the respective boundary and/ or initial conditions.
When the flow is viscous and laminar, a term appears in the equations of motion
that is related to shear stress which, according to Stokes' law of viscosity,
is
du
dy
where the subscript "l" denotes ''laminar, " y is normal to the direction of the
flow velocity, u, and ij, is the molecular viscosity (which is a property of the
fluid depending only on temperature under ordinary conditions). In compressible
laminar flow, when temperature variations are involved, the heat flux (following
Fourier's law of thermal condition) is shown as
^ dT
where k is the molecular thermal conductivity, which is also a property of the
fluid and usually depends only on temperature.
On the other hand, if the flow is turbulent, a kind of motion with random
irregular fluctuations is superimposed on the main flow. Since its fluctuations
are irregular, it is impossible to describe the motion in all details as a func
tion of time and space. However, one can study its average behavior using the
statistical methods pointed out by Hinze (1959). Therefore, in describing a
turbulent flow in mathematical terms, it is convenient to separate it into a mean
motion and a fluctuation, such as
u == u + u
T = T + T'
26
where the "bar" indicates mean quantities and the "prime" indicates the fluc
tuation from the mean.
The presence of random irregular fluctuations of eddy motion of the flow
appear as apparently increased viscosity or diffusivity (Defant I961) . In tur
bulent flow, therefore, additional Reynolds stresses must also be considered in
the general governing equations of the motion.
Analogous to the coefficients of molecular viscosity, Boussinesq (1877)
introduced a mixing coefficient, A^, for the Reynold stress in turbulence ex
pressed in terms of the mean velocity gradients in the flow field,
^ _ du
t " " ^ u^v' " t dy
where the subscript t denotes turbulent flow. In a similar definition to that
for molecular diffusivity, Schlichting (1968) postulated that the turbulent
heat flux is
t P q qy
The exchange coefficients of turbulent momentum and turbulent heat flux,
A^ and A^, both have the dimension of a viscosity, \i (g/cm/sec).
If the turbulence is isotropic (i.e., the relation between stresses and
strains is independent of axis) then A. corresponds to the molecular viscosity
\i, A corresponds to the molecular heat diffusivity (k/c ); and e = At/p, and
the eddy kinematic viscosity corresponds to molecular kinematic viscosity
V = [x/p* Turbulence shear stress is expressed as
du
T = p € — .
t dy
The eddy kinematic viscosity, e can be expected to occur as a constant only
if the turbulence field is homogenous. Generally speaking, the eddy viscosity
is 10^ to 10^ greater in magnitude than molecular viscosity, except in the vis
cous sublayer very near the solid boundary. The latter may be neglected in
many cases to a good approximation. As pointed out by Lamb (1932), in describing
the wind produced current, \i is replaced by A. .
By the techniques mentioned above, it is possible to attack geophysical
fluid problems which are usually turbulent in nature. In this way we can probably
predict or evaluate the timemean values of these complicated turbulences whose
complete theoretical formulations have so far proved impossible (Lamb I952).
27
Most of the geofluid modelings are treated in this way. For example^ previous
work by Sverdrup (19^7)^ Stommel (19^8)^ Munk (1950)^ Munk and Carrier (1950)^
Bryan (I965), and Birchfield (I967) all follow this kind of approach.
Nevertheless not all the turbulence is homogeneous and isotropic. The
eddy viscosity usually does have directional preference, Eddy viscosity be
comes effective only if there is some flow in the fluids and its value depends
on the velocity distributions and on the characteristic length of the flow field.
As shown in Boussinesq^s hypothesis and von Karman*s similarity theory
t dy 'dy
where L is the Prandtl's mixing length. Therefore A^ or e is not a property of
the fluid itself^ whereas the molecular viscosity^ \i^ is. It is still extremely
useful to introduce the eddy viscosity in relating the Reynolds stress to the
mean velocity gradients of the flow. Though A, (or g) is not a physical con
stant^ but varies from one part of the flow to another, estimation of the magni
tude of the average eddy viscosity is still the key factor to the success of
geophysical fluid modeling. In Lake Michigan, though the mean velocity is in
general less than l/2 m/s (Ayers^ et al. 1958)^, the flows are turbulent in
nature. Therefore, it is important to estimate the order of magnitude of the
eddy viscosity and the eddy diffusivity in order to facilitate geophysical fluid
modeling of Lake Michigan.
THE MAGNITUDE OF EDDY VISCOSITY AND EDDY DIFFUSIVITY IN LAKE MICHIGAN
In order to evaluate the magnitude of the eddy diffusivity and eddy vis
cosity associated with the thermal current structure of the late spring/early
summer circulation in the Lake Michigan model, a series of experiments were
carried out to verify order of magnitude values for the eddy diffusivity and
eddy viscosity to be used with further theoretical investigations.
An Estimate of Eddy Viscosity by Reynolds Stress and Vertical Velocity Gradients
In August 1967^ a pair of current meter buoy stations were set 8 miles off
shore from Sheboygan, Wis. (see p. 12 of this report for a description of the
buoy stations and their locations) »
During the period of the current measurement, the research vessel INLAND
SEAS was operating around the stations (during the daytime) whenever weather
permitted. Wind, as well as other weather data, was automatically recorded on
a digital meteorological tape recorder. The anemometer and vane that feed into
the recorder are on the mast of the ISLAND SEAS ik m above the water surface.
The Captain ^s log and hourly meteorological observations recorded wind from a
28
lower vane and anemometer that are about 11 m above the surface. Using the
analyzed recorded mean wind versus height plotted on a log scale with the mean
velocity zero at about 1 cm height (Floyd C. Elder, personal communication) the
10m level wind was interpolated with reference to the Captain *s log. The 10m
daily winds are shown in Table Jl. The vertical velocity gradients in the
direction of surface current were deduced from two current meters separated
vertically by 10 m. According to Deacon (I962), the Reynold stress is
= A ^
t dy •
— 3
Then using p^ = 1 g/cm^ and p^^ = I.25 x 10 g/cm^ the eddy viscosity can be
calculated. The data^ wind, wave conditions, absolute vertical velocity gra
dient, shear stress, and calculated eddy viscosity are all shown in Table 5I.
The values of Reynold stresses obtained here using Deacon's formula are in good
agreement with the mean value of stress by Elder and Soo (1967) at The Univer
sity of Michigan Research Tower on the east side of the lake.
An Estimate of the Upper Bound Value of Eddy Viscosity by the Result from
Statistic Theory
It is quite well known (e.g., Ayers, et al. I958) that a large clockwise
eddy exists in southern Lake Michigan between Grand Haven and Michigan City.
This large eddy has the length scale of a little less than half the dimension of
the width of the lake. Bowden (196^1) has shown the effective eddy viscosity is
proportional to the V5 power of the length scale, L, of the particular eddy
under consideration. In Lake Michigan, let this largest eddy have average
length scale of L = ^0 km or L = i+ x 10^ cm. Then, according to Richardson's
formula, the eddy viscosity will be of the order lo'^ cm^/sec which is the upper
bound of the numerical values of the eddy viscosity in Lake Michigan.
An Estimation of Vertical Eddy Viscosity from Current Measurement Near the Bottom
From 11 May to 11 October ±9^7 y a tripod supported pendulum current meter
was installed about 25O m offshore, Tl/2 miles south of Benton Harbor in a
depth of 4.5 m. The pendulum was suspended only 20.3 cm above the sandy bottom.
The monthly mean current is shown in Table 52. Using Lesser 's (195I) roughness
length for sand bottom, z^ = .I5 cm, the friction velocity near the bottom can
be determined through
z + z
\ = V[in(^)]"'
o
29
(D
H
•H
B
00
P
q5
VD
ON
■P
S,
(M
I
^
O
•H
■P
cd
■P
w
o
P
I
+:>
(L>
O
o
o
■P
•H
Ul
O
O
m
>
cd
t>.
o
rQ
CD
CO
>?
Id Ch
H O
H
I
pq
>.
+:>
•H
w ^>
o o
O 0)
w w
o
•H tsl \>.
(M
> <; CM^
VO
g
>= o
T^ v_.
^
w
H
Ul
*
O W
CD W
CD P
Ph a:!
CO >^
0.
O
•H
P
O
(D
CD
■P
cd
cvi
I
VD
CJ\
CM
\\ c^j
H
i
bO
O
H
CM
CO
CM
CM
H
00
H
CM
CM
CM
I
O
H
VD
H
CO
CM
's
1
O
H
X
XI
00
00
H
LTN
CM
CM
H
CQ
o
CM
t
H
CM
r^
CO
CM
CO
I
O
H
H
CM
CO
m
Ul
CQ
»
CO
(D
ir\ CD
hTN CD
g
LTN
CM
CD
>
1 !>
1 t>
H
1 1>
1
fe
cd
 cd
 Cd
cd
 Cd
o.
Cd
>
N^ ^
H ^
o
K^ >
H
^
o\
CO
bO
CM
X
(D
•
a
H
o
cd
>
"So
CD
CO
bD
cd
b
^
H
CD
1>
X!
cd
K^
CM
30
where u is the mean velocity, z is the height from the bottom during measure
ment, and Kq = O.Ul, The monthly mean friction velocity for z = 20.0 cm is
shown in Table 52. The neutral stability condition generally persists near
the bottom. Then the vertical eddy viscosity can be estimated as shown in
Table 52 (Bowden I96J4) .
TABLE 52. Estimated vertical eddy viscosity by current measure
ment near bottom.
Vertical eddy
Mean current Friction velocity viscosity
(cm/sec) (cm/sec) (Az = k u z cm^/sec)
Month
o^
May
59.5
June
51.7
July
1+7.2
Aug.
58.5
Sept.
56.1
Oct.
69.8
k.6k
2.krj
5.68
I+.56
U.58
5.UU
UO
20.7
50.9
38.5
56.0
1+5.7
Average value
33 cnf /sec
An Approximation of the Vertical Eddy Viscosity by Surface Wind
Wind data from both the two currentmeterbuoy operation off Sheboygan
and the wind and current station outside of the Benton Harbor beach are used
for approximating the vertical eddy viscosity in Lake Michigan. In the two
currentmeterbuoy operation off Sheboygan, wind data were measured aboard the
R/V inland seas. At the current meter station outside of the Benton Harbor
beach, the wind was measured Just above the water. From Sverdrup, et al. (I96I+),
depending on whether the magnitude of wind was greater or less than 6 m/s, the
vertical eddy viscosity was approximated as shown in Table 33. Notice that
the vertical eddy viscosity is smaller when the wind is averaged for a longer
period of time .
31
TABLE 33. Estimated vertical viscosity by wind speed.
Vertical eddy
Time Wind (m/s) viscosity, A^ cm^/sec Remarks
IT Aug.
1967
9.5
18 Aug.
1967
9.8
20 Aug.
1967
5.7
21 Aug.
1967
5.^
25 Aug.
1967
6.2
2k Aug.
1967
h.6
581
165
100
Wind averaged
over 6 hours
Average value
206 cm^/sec
May 1967
3.5
June 1967
2.5
July 1967
5.1
Aug. 1967
5.2
Sept. 1967
2.8
Oct. 1967
k.l
Nov. 1967
k.k
i;^
Monthly mean
wind
16
50
55
22
TO
87
Average value
^3 cm^/sec
A Measurement of the Eddy Diffusivity by the Dispersion of Marked Particles
In a series of experiments conducted in Lake Michigan, the fluorescent dye
rhodamine B was used as a tracer for diffusivity measurements. Rhodamine B has
peak ultraviolet fluorescence wavelength of 58O millimicrons. In May I968, a
fluorometer with only one sample cell was used. Because in this season the
phytoplankton was so abundant, and since they have an ultraviolet fluorescence
frequency Just a little lower than that of rhodamine B, the fluorescence of the
peak concentration of rhodamine B dye was not distinguishable from the background
after an hour or so. In all subsequent experiments, a comparative measuring
fluorometer with two cells was used. The rhodamine B mixed with natural lake
32
water flows through the dye sample cell, while the standard cell contains natural
lake water. The relative difference in fluorescence between the two cells is
measured in order to obtain reliable data in the possible presence of any inten
sity variation of ultraviolet light lamp,^ systematic changes in the characteris
tics of photocell, and differences in the amount of naturally occurring fluorescent
substances in the water (Noble and Ayers I961) • Both instantaneous source and
continuous source measurements have been made in I968.
Instantaneous Sources . Three runs of instantaneous sources of marked
particles for diffusivity studies were conducted. These were 1 and 8 miles off
Grand Haven in May, and 5 miles off Milwaukee in June. Rhodamine B dye was
diluted with methanol to have a specific gravity very close to the lake water.
Generally, one part of stock solution of rhodamine B is Uo^ acetic acid solu
tion was mixed with 51/^ parts of methanol. The dye mixture was gently poured
out on the surface lake water at the upwind side of the boat in order to allow
the boat to drift away with least interference to the dye patch. The first and
second runs were investigated aboard the r/v MYSIS and the third run aboard the
R/V INIAKD seas. During the first 20 minutes after the source had been instan
taneously releaseid, the diameters of the expanding dye patch were estimated by
eye with the ship's length as reference. At 20minute intervals, the boat was
used to coast through the dye patch longitudinally and laterally at the lowest
possible speed. For each transect, fixes were obtained from shore objects when
entering and leaving the patch, along with the time required to traverse through
the patch. Continuously running back and forth across the dye patch, the con
centration of the dye patch was continuously recorded from the fluorometer on an
AZAR recorder. The concentration distribution of a typical pass through the
patch is shown in Figure J"! Most generally, the dye concentration distribution
showed a form of Gaussian distribution. The diameter of the patch, longitudinally
or laterally, is properly approximated as k times the standard deviation of the
patch from its mean position. Batchelor (1955) shows that (assuming the fluid
is in a uniform, homogeneous, and steady, mean velocity field) at the initial
stage of diffusion process, the effective diffusivity (which is proportional to
the rate of change of the square of standard deviation with time) is proportional
to the first power of the elapsed time; in the intermediate phase, it grows as a
V5 power of the patch size; and at the final stage, as the standard deviation
grows as l/2 power of the elapsed time, the diffusivity should asymptotically
approach a constant. The general information and the value of effective diffu
sivity are shown in Table JU. The growth of dye patch versus time is shown in
Figures 52, and 55.
Continuous Source . Four runs of continuous (in time) point sources of marked
particles were conducted near the meteorological research tower 1 mile off the
east shore in Lake Michigan. Ten gallons of rhodamine B and methanol mixture
were contained in two 5gallon carboys connected together by plastic tubings as
shown in Figure 5^. The continuous point source was released from the apparatus
through a rubber tubing which had its outlet about 20 cm above the water surface.
The rate of discharge of the dye was controlled by a screw clamp on the tubing,, and
was usually allowed to flow at a rate of 5 gallons/hour, or approximately h g/sec.
In a homogeneous fluid with a steady uniform current, the continuous dye source
was observed as a slender plume extending downstream.
55
r «' 4'
Distance (in unit scale)
6^UT
Figure 51 A typical distribution of dye concentration across a dye patch.
TABLE 5i^■. General information and effective eddy diffusivity ©f dye
patch studies.
Run
Run 1
Run 2
Run 5
Date
1 May 1968
2 May I968
27 June 1968
Location
1 mile offshore
8 miles offshore
k miles offshore
near Grand Haven
near Grand Haven
near Milwaukee
Depth
20 ft
ko ft
65 ft
Duration of run
1150  1?10 EST
1500  IU50 EST
0910  12if0 EST
Current velocity
cm/sec
11.6
5.5
7.6
Wind
Direction
m
SW
SSW
Speed (kts)
15
1520
5
Lake condition
25 ft waves
1 ft wave
2 ft waves
Effective eddy dif
fusivity (cm^/sec)
Lateral
Not available
lli^O
7^5
Longitudinal
Kot available
^131
^019
Remark
Very poor visi
bility
34
l/>00<
600
300
too
60
30
,T^ .rife
f
/
y
/
/
/
/
/
(a)
T6 '1.2 • ■3.0' ' ■ ' "6
Time (minutes)
1,000
600
30 0J
I
o
too
60
/
y/z
/
/
/
/
/
X^
/
/
X/
(b)
n:! ^ ' *3^' ' ' ' '6 Hi ' ' '30 60 ^20"^' '
Time (minutes)
Figure 32. Lateral growth of the dye patch with time,
(a) 2 May I968. (b) 27 June 1968
55
I,000
/*
600
/
300
/
/
/
4J
0)
/
w
/
^
X
100
/
/
/
SO
SO
10
(a)
600
100
Time (minutci}
/
X
/
/
X
/
/
/
/
/
/
/
/
x/
/
/
(b)
. Ti—" — ■ "sd ■ ■ ' V. ■ "« ■ 'so' ' ■ ' "60 "So"" — '
Time (minutes)
Figure 55, Longitudinal growth of the dye patch with time,
(a) 2 May 1968. (b) 2J June I968.
56
Air Intake
Screw Clapper*
Figiire 5^. Continuous source apparatus.
The direction of the dye plume was reported by the observer using a theodo
lite mounted on the meteorological research tower 57 ft above the water surface.
One or two buoys were dropped as markers on the edges of the dye plume at right
angles to the direction of flow. The appropriate distance of the marker from
the source depended on the mean velocity of the current which was estimated from
the drogue study set earlier.
The r/v MYSIS crossed the dye plume back and forth at 10minute intervals
at the lowest speed she could make (about k.^ mile/hour). As soon as the dye
plume passed through the markers, samples of water were taken through a constant
flow fluorometer that fed into an AZAE recorder. Communication between the
tower and the boat was maintained through a VHF transceiver. Whenever the boat
entered and left the dye plume a "MARK" signal was transmitted to the theodolite
observer. The observer obtained the azimuths of both edges of the plume and the
elevation from the tower to the boat.
The horizontal distance between the boat and the tower was obtained from
the elevation angle as well as by fixes from land objects. The angle subtended
by the edges of the plume is the difference of the two azimuths. The arc length
of the plume at this distance from the tower was arbitrarily taken as four times
the lateral standard deviation unit. At each return crossing, the vertical
diffusion was measured by lowering the fluorometer intake at half meter inter
vals until there was no indication of the dye. Half of the maximum vertical
displacement was taken as the vertical standard deviation unit.
57
At the last few minutes before the exhaust of the point source^ continuous
crossings were conducted at various distances along the plume from the tower
to approximately II/2 miles downstream. These crossings measured the asympto
tically steady diffusion of the plume.
These continuous point source diffusion studies were carried out on 50 July
and 12 August I967. On 50 July, sampling was discontinued about 2 hours after
the experiment was initiated, because of high seas. The dye plume was contin
uously sampled for 5 hours on 1 August, and during the sampling period the cur
rent shifted its direction of flow about I5 degrees. Two sampling runs were
conducted on 2 August under a very steady current condition. The general infor
mation for each run and the value of diffusivity are shown in Table 55. The
lateral growth of dye plume versus time at a constant distance from source is
shown in Figure 55. The lateral growth versus distance under asymptotically
steady condition is shown in Figure 56.
A Measurement of the Horizontal Diffusion by Drogue Studies
A drogue is made of k pieces of light sheet metal (e.g., aluminum alloy)
approximately 1 x 2 m in dimension, linked together in the shape of a cross
and suspended by chains in the water by a float which supports a radar antenna
above the surface (similar to Farlow I965). It is frequently used for measuring
the mean current of the flowing layer at the drogue depth. The equivalent radius
of the drogue is about 1 m and, based on its geometry and the corresponding
Reynold numbers (10^ to 10^ in the Great Lakes), its time constant is estimated
to be 100 sec.(Okubo and Farlow I967) . Though the smallscale turbulence may
be averaged out, the drogue is adequate for measurement of largescale turbulent
diffusion for eddies larger than 10 m in scale which represent the energy
containing eddy scales in the Great Lakes (Csanady I965; Okubo and Farlow 1967) .
On 27 June I968 a hexagonal pattern of surface drogues, with the center
one very close to a fixed buoy, was set about 5 miles offshore near Milwaukee
in Lake Michigan. The original distance between each pair of drogues was l/2
mile. After setting the pattern, the r/v INLAND SEAS continuously fixed each
drogue for the following 6I/2 hours. Similar to the result mentioned by
Okubo and Farlow (1967), the growth of the pattern size with time was not
clearly noticed. The diffusion constant was deduced from the second moments
about the mean of the haxagonal pattern. At any particular time, the instan
taneous positions of all drogues were interpolated from their respective pre
vious fixes. The general information on the experimental run and the value of
eddy diffusivity are shown in Table 56.
COMPARISON OF THE VALUES OF EDDY VISCOSITY AND EDDY DIFFUSIVITY WITH OTHER
STUDIES
In Lake Michigan, there is little data about eddy viscosity or eddy dif
fusivity. Some figures from the ocean or other lakes are available as reference «
58
CD
O
U
g
w
§
•rl
P
o
o
O
>:.
P
•H
t>
•H
W
;:^
Ch
Ch
•H
<D
CD
>
P
O
CD
Ch
Ch
CU
cd
s^
O
o •
Ch >5
P
m
cd
0)
es ft
pq
CD
•s
EH
iH
CD
CO
^
Sh
H
CX)
^
•H
i)
MD
Ch
•H
d
ON
Ch
cd
•n
^
LTN
H
§
§
•H
H
•
CD
^
bO S
1
d
tiD
H
P
•H
?H
r^
CD
Ch
j;^
<
S
C^J
W
M
H
J2J
CD
5
cd
CM
«
CM
H
CD
CD
S
h^
\s\
H
EH
M
fH
pq
XI
c6
a
Cd
CX)
w
Hi
•H
W
LTN
VD
Ch
•H
^
C^
Ch
cd
•\
^
On
H
g
g
»
CD
X
bD 1^
1
hr\
tiD
H
CD
P
^
'd
U
M
CD
Ch
5;;^
<
n
cd
OD
^
jj
CD
5
Cd
CX)
p:^
CVJ
H
CD
a?
S
H^
Lr\
EH
CO
X
S
W
^
cd
a
cd
00
Ul
h:)
•H
bD
1^
VD
Ch
•H
hf^
o^
Ch
cd
•\
^
hCN
H
Cd
H
•
S
bD :^
1
CM
bO
H
CD
P
jii
•H
^
r^
<D
Ch
\r\
a
«^
S
cd
w
^
H
j:^
CD
5
•^
pc;
H
H
CD
s
h^
ir\
EH
CO
>:!
a
H
00
^
cd
a
cd
VD
w
iq
•H
bO
ON
Ch
•H
CM
H
Ch
cd
•\
^
CM
0.
a
H
!>5
cd
•H
H
CD
S
bO S
1
H
H
CD
P>
h)
^
M
CD
Ch
f^
^
cd
CO
cM
a
jj
CD
jj
Cd
r—
P^
rr\
H
0.
I^
i:
Lf^
r—
a
•H
P
CD
Cd
P>
Cd
P
iq
■p
Pi
CD
Ch
o
o
p
cd
Pi
00
CM
CO
CM
00
CM
l:^ 00
^ 00
p
o ^
O TJ
1:1
CD
•H
P
CD
CD
ft
CO
CD
>
Cd
■p
Ch
CM
I
cd
p
Ch
I
CM
CO
CD
!>
cd
15
p
Ch
\r\
w
1
CO H
hr\
o
•rl
P
•H
t;}
a
o
o
CD
cd
P
CM
CO
(D
!>
cd
^
p
Ch
H
CM
I>
1
hCN
rTN
CO
Lr\
X^
CD
bD
Q)
g
;:s
H
•H
rQ
P
0.
cd
CD
H
bD
•H
0.
bD
Cd
•H
a
%
CO
H
P
cd
P
CD
^
X
s::!
CD
P>
CD H ^N
> > O
•H 'H
p> CO CO
O p( \.
CD Ch c\j_,
Ch Ch g
Ch H o
H Id ^^
CD
p:;
39
1.000
900
800
700
600.
500
u
(U
41
400
>>
300
200
100
/
/
/
/
60
120
(a)
'240 ' ^363 ' ■ '600
400
300
200
100
90
801
70
60^
50
18
/
^
.X'
'30 ■ '42' • '60
Time (minutes)
(b)
'120 ^I80"
Figure 55. Lateral growth of dye plume at a constant distance
from the source, (a) 1 August I968. (b) 2 August I968.
ko
\
\
O v
CD
O
1 I r
(53aj) aouBiSTQ pasnjjTQ T"35B1
\,
\.
\.
is
\.
\.
\
\
\
1 i i § i i"
« rt N M  =
(53aj) 33UB5STa pasnjJTQ TW35B1
o
o
o
cd
P
m
P
•H
^ •
B 00
3 vo
H ON
ft H
*i
o <ci
^ CM
■P
o ^
H ^^
cd
0) CO
+3 vo
cd ON
h^ H
P
• CQ
I
<t)
^ H
a
§,.
•H (d
s s § s
♦ •" lO «•
(3333) souBqsja P3snjjTa TBa35Bi
41
TABLE 56. General information and eddy diffusivity of
drogues study in Lake Michigan.
Date:
Location:
Depth of water:
Depth of drogues:
Duration of run:
Current velocity:
Wind:
Direction
Speed
Lake condition:
Mean initial separation:
Final deviation from mean:
Effective eddy viscosity:
27 June 1968
5 miles off Milwaukee
65 ft
5 ft
I5OT  1850 EST
7.6 cm/sec
N
6 kts
1/2
1 ft
waves
2578
ft
3UOO
ft
5A X
: 10*
cm^/sec
In modeling the Gulf Stream, Munk (1950) estimated the lateral eddy viscosity
of 5 X 10'^ cm^/sec. Stommel (1955) estimated the eddy viscosity of Gulf Stream
near the Florida Straits to be 2 x 10^ cm^/sec. Defant (I96I) determined the
exchange coefficient of lateral largescale turbulence in ocean phenomena as
5 x 10'^ cm^/sec.
In general, the vertical turbulence scale is much smaller in dimensions
than that of the horizontal. Therefore the magnitude of vertical eddy viscosity
is much smaller. According to Defant (1961), Sverdrup obtained the value of
the vertical eddy viscosity of 75260 cm^/sec in the North Siberian Shelf re
gion; Fjeldstad estimated the vertical eddy viscosity in the Caspian Sea of
022if cm^/sec, and Suda measured the value of 6807500 cm^/sec for the Kuroshio
and of 1501^^60 cm^/sec in the Japan Sea.
For vertical eddy diffusivity, some typical values have been summarized by
Defant (I961) as: Mediterranean, k2 cm^/sec; California current, 50^0 cm^/sec;
Caspian Sea, 2l6 cm^/sec; Equatorial Atlantic Ocean, 520 cm^/sec. Csanady
(1965,196^,1966) has measured the eddy diffusivity in Lake Huron and Lake Erie,
and found horizontal diffusivity of the order of 5OO25OO cm^/sec in Lake Huron
k2
and, 10005000 cm^/sec in Lake Erie. These values are much greater than verti
cal eddy diffusivity, which is only of the order 0.110 cm^/sec. Okubo and
Farlow (1967) measured the horizontal diffusivity for large eddies in Lake
Michigan and Lake Erie using drogues and found the eddy diffusivity of 5 x 10
to 6 X lO''^ cm^/sec.
According to our estimation, the vertical eddy viscosity is on the order of
55 to 2 X 10^ cm^/sec in Lake Michigan. If the vertical eddy viscosity is
usually two orders smaller than the horizontal value (Defant 196l)jj then the
horizontal eddy viscosity should be in the range of 3.5 x 10^ to 2 x 10^ cm^/sec,
which is within the upper limit of lO"^ cm^/sec as we have estimated. The
lateral eddy diffusivity we found using marked particles is 7501200 cm^/sec for
small eddies and ^.k x 10^ cm^/sec for larger eddies. The effective longitudinal
eddy diffusivity for small eddies is in the order of ^ x 10^ cm^/sec. All the
data we used appear quite consistent and in good agreement with that of others.
CONCLUSION
In modelling a geofluid problem^ the eddy viscosity and the eddy diffu
sivity are of critical importance in order to predict the natural current or
wave phenomena with a similarity to laminar flow. In Lake Michigan, though the
mean current velocity is small in general, the flow field is turbulent in
nature. In solving or explaining the flow pattern of the mean lake current, it
is possible to use the governing equations of laminar flow with eddy viscosity
and eddy diffusivity in place of the molecular viscosity and the molecular dif
fusivity* A series of experiments have been conducted in order to estimate
these two values. The data presented are quite consistent and in good agree
ment with the data reported by other investigators. The vertical eddy viscosity
in Lake Michigan is in the range of 1 to 10^ with a mean 10 cm^/sec. The
horizontal viscosity is in the range of 10^ to 10"^ with a mean value of 10^
cm^/sec. The eddy diffusivity may reach the same magnitude as the viscosity
but it is in general smaller. A typical mean value for vertical eddy diffu
sivity is 5 cm^/sec^ and for the horizontal eddy diffusivity 10^ cm^/sec.
Though few experiments are far from enough to state the characteristic values
of turbulent phenomena, the data presented together with values reported by
others can be used as good reference values for Lake Michigan.
h3
References
AYERS, J. C, D. C. CHAKDLER, G. H. LAUFF, and C. F. POWERS. I958. Currents
and water masses of Lake Michigan. University of Michigan^ Great Lakes Res,
Div. Pub. No. 5, 169 p.
BATCHELOR^. G. K. 1955 The theory of homogeneous turbulence. Cambridge Uni
versity Press, Cambridge, Mass. 197 P*
BIRCHFIELD, G. E. 1967. Winddriven currents in a long rotating channel.
Tellus 19(2):2i52i9.
BOUSSINESQ, J. 1877. Theorie de l^ecoulement tourbillant. Mem. Prec. Acad.
Sci., XXII, I46, Paris.
BOWDEN, K. F. 1964. Turbulence. P. II50 in Oceanogr. Mar. Biol. Ann. Rev.,
Vol. 2. H. Barnes, ed. George Allen and Unwin., Ltd., London.
BYRAN, K. 1965. A numerical investigation of a nonlinear model of winddriven
ocean. J. Atmos. Sci., 20:59^606.
CSANADY, G. T. I963. Turbulent diffusion on Lake Huron. Jour. Fluid Mech.,
17:560585.
. 1966, Accelerated diffusion in the skewed shear flow of lake
currents. J. Geophys. Res., 7I: 4llJ^20.
DEACON, E. L. 1962. Aerodynamic roughness of the sea. Jour, of Geographic
Res., 67:51575172.
DEFANT, A. I961. Physical oceanography. The MacMillan Co., New York, Vol. 1,
729 p.
ELDER, F. and H. K. SOO. I967. An investigation of atmospheric turbulent
transfer processes over water. University of Michigan, Great Lakes Res. Div.
Spec. Rep. No. 29, k6 p.
FARLOW, J. S. 1965. A field technique used for horizontal diffusion studies
in Lake Michigan and Lake Erie. University of Michigan, Great Lakes Res* Div.
Pub. No. 15, pp. 299505.
HINZE, J. L. 1959. Turbulence. McGrawHill Book Co., New York, 586 p.
LAMB, H. 1952. Hydrodynamics. 6th ed., Dover Pubs., New York, 758 p.
LESSER, R. M. 1951* Some observations of the velocity profile near the sea
floor. Trans. Am. Geographic Un., Vol. 52, pp. 207211.
kk
MIMK, W, H. 1950. On the winddriven ocean circulation, J. Meteor,^ 7(2):
7995
and G. F. CARRIER. 1950. The wind driven circulation in ocean basins
of various shapes. Tellus 2(5) : I58I67.
NOBLE, V. E. and J. C. AYERS. I96I. A portable photocell fluorometer for dilu
tion measurements in natural waters. Limnol. Oceanog. 6:U57^6l*
OKUBO, A. and J. S. FARLOW. 19^7. Analysis of some Great Lakes drogue studies.
Proc. lOth Conf. on Great Lakes Res., Internat. Assoc, for Great Lakes Res.,
pp. 299508.
SCHLICHTING, H. I968. Boundary layer theory. 6th ed. Translated by J. Kestin.
McGrawHill Book Co., New York, jkj p.
STOMMEL, H. 19^8. The westward intensification of winddriven ocean currents.
Trans. Am. Geophys. Union 29(2) : 202206,
SVERDRUP, H. U. 19^7. Winddriven currents in a baroclinic ocean; with appli
cation to the equatorial currents of the eastern Pacific. Proc. Nat. Acad.
Sci. 55(11): 518526.
^5
IV. MJMERICAL EVALUATION OF STERN ^S CIRCULATION MODEL
James H. Saylor
DEVELOPMENT OF STERN^S MODEL
In 1965^ Stern presented a model for the interaction of a uniform wind
stress with a geostrophic vortex. The results of these calculations showed
that a geostrophic vortex subjected to a uniform wind stress would move, as an
entity, with an Ekman drift, provided that the thermocline was "soft" and that
internal waves could be generated at the center of the vortex. Stern* s results
have been evaluated in terms of wind stresses and vortex dimensions that might
be expected to occur in Lake Michigan (Noble 196?). The results predict a very
broad range of internal wave periods that are within the limits of values
measured for Lake Michigan (Mortimer 1965).
The starting point for the theoretical model is the assumption of a hydro
static, twolayer model; upper layer of density p and lower layer of density p^;
no friction; and the fluid incompressible.
A
A
Dv J ^ VP
A J1
—  + fkxv+ — = •
B'^*^
Dt p
(1)
where H = depth of upper layer.
and
p  p
o
D d ^
^ A A A A A
U = V + Wk = Ui + vj + Wk
A
VU =
Let f/h = (S/Sz) G where 9 is permitted to be a lateral wind stress only. ^In par
^ / \ A ^
ticular, take (x,y,z = o,t) equal to the surface wind stress, t. Let U =
Uq + Uj^, where 11^ is taken as the geostrophic component of (l).
k6
A ^ h
U = V + KW
o o o
A / A
u = V + m^
b b b
Thus
A
D V
o o
Dt 'op
A /I VP /
+ fk X V + — =  g'k (2)
o
where
A
VU =0 ,
o
A A
and whence^ Uj^ is readily defined as the component of U due to the surface wind
stress, T. Call this part the "Ekman drift." Boundary conditions are established
by requiring no displacement of the free surface^ i,e,^ ^'U^ =  k^^U^ =
W^(x^y^z = O)^ and by requiring at the fluid interface that the "Ekman drift" is
limited to the upper layers^ and that the wind current joins the geostrophic
component smoothly Just above the interface,
^A
U^(x,y,z = H,t) =^(x,y,z = H,t) = ,
which implies
5. ..^ ^ ^^ n
9(x,y,z = H,t) = ^ = ^ =
A/ s A
9(x,y,z = 0,t) = T
Subtracting the geostrophic component, equation (2), from the momentum equation,
equation (l), gives:
A
which can be written as:
b^AAAAA Aa^G ,,
T— + V • Vv + V Vv + V Vv + f k x V = :r . ( 5 M
at boobbb bdz ^^ ^
hi
Taking the vertical component of the curl of equation (5')^ denoting:
if s k.V X U = kV X V
^o o o
A 4 4 A
Stern obtains,
where:
^ A /I /I A
2^ = kV X U = k*V X V
^b b b
(f + <f )^— + ^(kV X 0) = v^ViT + Z T. (4)
o' az 3z' b ^o i ^ ^
1=1
dz dx o b dz dry' o b
o D o b
Ts
Sz 5x dz ^ '
a
<^
^ .. ^f
Ts = (W + W^)^— + W^ .
o b dz b dz
Integrating equation (h) between z = H and the free surface:
/^ (f + ^ )^ dz + /^ (kV X 0)dz = /^ v^.V^'dz + /^ Z T.dz (5)
H ^ ^o^Sz H Sz^ ^ H b o H . , i ^^^
1=1
which becomes:
•(f + W )W^  kV X T + f V V? dz =  / Z T.dz . (6)
o H b o H . , 1
1=1
i+8
stern compares the magnitude of the terms in equation (6) and shows that if the
Rossby numbers based on both the geostrophic and "Ekman" currents are small,
i.e. y
D ^ O _
the / T. terms are of negligible magnitude. Defining the winddriven volume
transport as;
/^ v^dz = M =   k X T
(which is the Ekman transport in the conventional sense).
(f + ^ )W^  kV X T + MV^ = (6')
where A a h
^ ^^ k X TV6
To include the effect of f varying with latitude, use the betaplane approxi
mation, which simply replaces the V<^ term with V(f + ^ ). So, the general so
lution is written:
which may be expressed as:
«* = 'ff^ (T")
where
For the vorticity of the geostrophic component in the upper layer we have the
relation:
^9
Integrating this expression between H and the surface:
,0d,. .,, ,0 ^
Ih H(' ^ i' = /^ S^  '
A A
dt
which may be written as:
d(f + ^ )
dt H
, f rr ^ X T f dH
f + £' ~ H dt
(9)
Equation (9) is the generalized equation of motion for the geostrophic component
in the mixed "Ekman layer."
For a constant, uniform wind stress (v x r = 0) and for a limited area in
which f ^ constant, equation (?') reduces to:
(k X T)*Vf
^^=77— F^ (10)
^^o
o
and the kinematic equation, equation (9)i becomes:
H dt " H ( f + f' )^
dt
. (11)
Stern's analysis with the necessary conditions and restrictions, describes the
interaction of a localized, geostrophicallybalanced flow (i.e., within the
confines of the simple, twolayer model envisaged) with a winddriven "Ekman"
current. The vertical velocity (W^) generated by this interaction is propor
tional to the gradient of the geostrophic vorticity (V^), a quantity which is
conserved in the upper layer, and to the wind stress, t. In this manner, it
is seen how the vertical velocity transfers momentum across the interface to
the lower layer, and how there is a tendency for the vorticity of the water
column to increase due to the advection of neighboring vorticity by the mean
Ekman drift.
50
APPLICATION OF STERN'S MODEL TO LAKE MICHIGAN
If the wind stress (t) and the coriolis parameter (f) are independent of
position, the last term of equation (ll) reduces to:
^(kXT)V^^ ^k.(V^^XT) k.Vx(^^fT)
a k.v X (— )
which is equivalent to a force per unit mass of  ^^r/Ef acting on the upper
layer.
Consider a twolayer, model (Fig. kl) with Ekman drift at right angles to
T and directed along the x  axis, i.e..
T = positive constant
y
T = 0.
X
^ ^_^ ^ ^ ^ ^ ^ ^ ^ ^  • ^ ^ X^
Vorticity <, I E^ /> .y \ + \
/// ^ /////// //> ^ '
Vorticity li a. ^ H^ '• ^ ' 5 H,  H,
Figure i.l. Schematic two layer model.
Let the positive displacement of the interface at z = Hi be small and given
"by Ti(x,y, t). The vorticity equation for the lower layer is:
d4 _ ^ ^
dt az •
Integrating between z = H^  Hi and z = Hi + r
,Hi+Ti d^ _ pHi+n ^ ^ ^^
iHiH2 dt ^^ iHiHe ^ hz ^^ '
dt Ha' 'z = Hi+Ti ,
51
^ ^z = Hi + n dt
Integration gives:
^*# = ^*it" '
t=0 t=0
^^(t) = Ul^llll . ^^(^^y^t=0) . (12)
The vorticity equation for the upper layer is:
'i^^i'^^^^t'"
where U = r/fHi = mean "Ekman" drift.
(15)
The interaction of an Ekman drift current with a geostrophic vortex can be
described in terms of a vertical suction velocity W^, which is proportional to
the gradient of the geostrophic vorticity and to the wind stress. This interac
tion leads to an increase in the total vorticity of the water column with time
because of the transfer of momentum across the interface. It was shown that if
T and f were constant, the result of the interaction can be described as advec
tion of vorticity at the rate of mean Ekman drift. Equation (I5) expresses this
fact. If the geostrophic vortex is centered at x = y = at time t = 0, at
which time an Ekman current is established, equation (I5) states that part of
the vorticity will stay at x = y = (represented by the interface disturbance),
while part is advected with the drift current.
An expression for the vorticity difference between the upper and lower
layers comes from Mar gules * formula
*ri =<&  f^vit) (Ik)
where
r^ ^ ^
^ ^""W
From (12), (15) and (lU):
/l. + ui)(^ .^Ivin) + — ^ = U — ^
^ ' "^ hi'^YLz f '^"'^ ' Hi ^ "" H2 ax ' ^"""^^
52
Separating equation (l4) into steadystate and time dependent components by
setting
Tl(x,y,t) = Ti^(x^y) + 7i^(x,y,t) ,
for the steady state component:
for the time dependent component:
Assuming a solution of the form:
ik(xCt)+i%
''t " °
(C = phase speed) then
\ = ^
vin^ = (k + i^)!^^
equation (I7) becomes:
^^^Hi_^^ + g' p(k^ + «^)] = U[l H g' ^(k^ + i^)] , (17.)
or,
§=[l+(lk^0"" (18)
where
K = (^)'/^ k
^, = (^)V^ i
55
4 1
Take f = 10 sec (an approximate value for Lake Michigan).
If k^ » Hg/Hi, C 21 U and short waves are advected at the speed of the
mean Ekman drift.
If k^ = i^ is small the phase speed is:
(1 + k^ + 4)Hi
and a surface vortex of radius ^(g* 3)^/^ will leave the origin in time:
^"2^ f^ ^ ' LuHi(l + 4 + 4) J •
To estimate the eddy diameters, ratio of propagation, and internal wave
periods involved, use the following parameters as indicated in Figure i— 2.
I  °
a
j^ = 20 m ^ = 0.9982 (20° C)
H2  100 m /tf© s^+A^= 0.9998 (8" C)
Figure k2. Specific parameters of the two layer model.
Let:
^:^ = ^:^ ^ f = 10" sec"
:' = (— )g « 1.57 X 10"^ m/sec^
T = 2 X 10' (m/sec)^
Then,
*T*
U = —  = 10 cm/sec
rHi
%
^ ~ ^[{1 H k^ + 4)Hi + hJ
eddy diameter = — ( ^g^ )"'"/^
TABLE i+l. Calculated:eddy diameters, phase velocities,
and interval wave periods of upper layer
in Lake Michigan.
K = h
Eddy
diameter,
km
cm/sec
Phase speed
c,
cm/ sec
Wave
period,
days
1
59
10
5.75
12.1
2
1.95
10
5.0
^.57
5
15
10
6.67
2.25
it
7.8
10
9.1
0.99
5
^.8t
10
9.55
0.59
Effect of Density Difference, Ap
Keeping other parameters fixed, let Ap = l/2 of former value, i.e., take
g' = 0.8 X 10 (cm/sec^). Considering the density change, the result is cal
culated as shown in Table 42.
TABLE 42. Calculated eddy diameters, phase velo
cities, and internal wave periods for
upper layer in Lake Michigan.
Eddy
k = i diameter, U,
km cm/ sec
Phase speed
cm/sec
Wave
period,
days
1
28.1
10
3.75
8.7
5
9A
10
6.67
1.63
8
5.5
10
9.55
O.J+5
55
Effect of T
The wind stress determines the speed of the Ekman drift. Within the limits
of (v[^fL « 1 then^ the only effect of t is to determine the rate of advection*
If the parameters of the basic model remain fixed, and t is reduced to l/2 of
its former value (1 x 10 (m/sec)^), the solutions are modified as shown in
Table k3.
TABLE ^45. Calculated eddy diameters, phase velo
cities, and internal wave periods for
upper layer with wind stress changed.
Eddy
^^ ~ ^^ diameter,
km
Phase speed Wave
U, C, period,
cm/sec cm/sec days
1
59
5
1.87
24.2
5
15
5
5.55
i^.5
8
k.Qj
5
U.76
1.18
Effect of Depth of Thermocline
Assume the same parameters and overall depth as basic model, but let Hi
10 m, H2 = 110 m, then the results will be as shown in Table hk.
TABLE kk. Calculated eddy diameters, phase velo
cities, and internal wave periods in Lake
Michigan with the depth of thermocline
changed.
K = h
Eddy
diameter,
km
cm/sec
Phase speed
c,
cm/sec
Wave
period,
days
1
5
i^5.4
14. U
20
20 1
4.28
12.6
11. T
1.52
For Hi + H2 = constant, the effect of decreasing Hi is to increase the
mean Ekman transport velocity (u), to increase the wavelength associated with
56
each wave number^ and to decrease the wave period (because of increased advec
tion rate) .
The range of values for internal wave periods calculated from this model
are within those reported by Mortimer (I965) and also within the range of perio
dicities of current speed fluctuations calculated by Verber (1965). However,
because of the limitations imposed by the assumptions of a simplified model,
these values cannot be considered to be definitive.
Stern* s model does provide mechanisms for the generation of phenomena
observed in field measurements, such as rotary currents from buoy stations,
rotary drogue motions, internal waves, and a transfer of wind energy into the
water circulation in a way that might explain the buoy station measurements re
ported in the previous section of this report.
It will be necessary to design future field experiments relating the three
dimensional temperature field to the threedimensional current field as a func
tion of time. These experiments will be timeconstiming and expensive, but the
numerical results presented in this paper indicate that there is a strong reason
to make definitive tests of the Stern model.
57
References
MORTIMER^ C. H. I965. Spectra of long surface waves and tides in Lake
Michigan and at Green Bay^ Wisconsin. ProCo 8th Confo on Great Lakes ReSo ;
University of Michigan^ Great Lakes Res. DiVo Pub, No, 15:50^1525.
NOBLE^ V. E. 19670 Evidences of geostrophically defined circulation in Lake
Michigan. Proc. 10th Conf. Great Lakes Res«^ Internat. Assoc, for Great Lakes
Res., p. 289298.
STERN, MELVIN E. 19650 Interaction of a uniform wind stress with a geo«
stropic vortex. DeepSea Research, 12:555""67o
VERBER, JAMES L. 1965» Current profiles to depth in Lake Michigan. Proc, 8th
Conf. on Great Lakes Res. ; University of Michigan, Great Lakes Res. Div. Pub.
Noo 15:56^571.
58
V THERMAL CURRENT STRUCTURE IN LAKE MICHIGAIT, ASSOCIATED WITH
THE SPRING WARMUP SEASON— A THEORETICAL STUDY
Joseph C. Huang
INTRODUCTION
The early spring condition of the lake lends itself most readily to mathe
matical analysis. During the winter^ the water mass of the lake becomes nearly
isothermal^ and has a threedimensional circulation pattern reaching to the bot
tom of the lake basin. In May^ "while the main water mass is still isothermal
and at its winter temperature of ^'^C or less^ the edge of the lake begins to
warm under the influence of a combination of factors, including warmwater run
off from the adjacent land. The nearshore boundary water warms to a temperature
about S^'C (and is isothermal to the bottom) in water depth of 20 m, and extends
a mile or more offshore. The edge between the boundary water and the main water
mass of the lake is characterized by a sharp temperature gradient (the "thermal
bar"), water color differences, marked differences in phytoplankton populations,
and strong geostrophic currents along the temperature gradient. As the season
progresses, the boundary water edge becomes more complex and often shows multiple
bands of water color and temperature 2 to ^ miles or even more offshore, parallel
ing the perimeter of the basin. The period of 28 May to 5 June 1967^ when this
process operated, has been described for Lake Michigan by Stoermer (I968) and
by Noble and Anderson (1968). This same process has been well documented for
Lake Ontario and Lake Huron by Rodgers (1965,1966).
In the Great Lakes area, the wind system is very much modified by the con
tinuous presence of a high pressure in the spring summer period and a low pres
sure in the fall winter period due to the meteorological effect of this world ^s
greatest freshwater masso The wind observed in the lake is quite different from
that of nearby land stations o Both the magnitude and the direction of wind are
frequently changing, seldom consistent over a long period of time^ The time
averaged net energy put into the lake by wind stress may not appear as a simple
windcurrent relationship. This is especially obvious during the lake warmup
season. In the spring, when the lake water is still relatively much colder than
the already warmed adjacent land, the air mass above the lake is frequently in
a stable condition. Wiresonde cross sections in Lake Michigan have demonstrated
the intense springtime inversion prevalent before the warming of the water sur
face (Bellaire I965) . At the Wisconsin shore, it has been observed that water
remains smooth even beneath offshore winds of greater than 10 knots (at 50 ft).
However, as far as the thermal energy connected with temperature difference be
tween land and lake and air and water are concerned, temperature plays a much
more constant role than the varying wind all the year around. During the lake
warmup season, temperature gradients, either horizontal or vertical, show per
sistence, which, in turn, suggests that the thermal energy input may dominate
the lake circulation dynamics,
59
Noble (1967) has indicated that the circulation dynamics of Lake Michigan
are probably dominated by geostrophic forces rather than the wind stress during
the summer season. It is the purpose of the present study to model Lake Michi
gan, considering the geostrophic thermal gradient force as the prime source for
the lake current generation, and to investigate the thermal circulation pattern
of the lake as a comparison with the general current reported by others (e.g.,
Ayers et al. 1958). The study has not reached the final stage yet. The model
may be refined by some modification much closer to the real simulation of the
lake. The preliminary result shows that the thermal circulation pattern in
Lake Michigan is in very good agreement with general current data reported
earlier (Ayers et al. 1958; Rodgers 1965) as well as our own investigations.
MODELING CONSIDERATIONS AND ASSUMPTIONS
In the geophysical fluid model of large scale motion, Coriolis effect must
be taken into consideration. During the spring warming season, a sharp hori
zontal temperature gradient always exists within a few miles from shore around
the perimeter of the lake. In a stable atmosphere, when winds are constantly
shifting above the lake, the thermal gradient force comes to equilibrium with
the Coriolis force. Generally speaking, only the geostrophic thermal gradient
system is apparently in operation during the spring warmup season. Equilibrium
between the temperature gradient force and the Coriolis acceleration produces a
component of horizontal velocity perpendicular to the direction of the imposed
thermal gradient in a rotating system. This is the geostrophic thermal gradient
relationship which is considered dominating.
In Lake Michigan, the length (about 500 km) is much greater than the width
(mean about 120 km). The wide Straits of Mackinac connect the northern end of
Lake Michigan to an even larger lake — Lake Huron. In the Lake Michigan model,
to the lower order of approximation, the flow field, away from the north and
the south ends, can be considered as zonally independent.
The Coriolis force varies with geographical latitude. In Lake Michigan,
its magnitude ranges from 0.^88 x 10 rad/sec in the southern end to 0.525 x
10 rad/sec in the northern end. The variation from south end to north end is
small. The Beta effect in the Lake Michigan model can be neglected.
An eastwest cross section of the central part of Lake Michigan is shown
in Figure 51. The deepest point in Lake Michigan is 28I m in the northern
half of the lake. The southern portion has a maximum depth of about I50 m.
In the Lake Michigan model, with a width of 120 km and depth of I50 m, a flat
bottomed trapezoidal crosssection as shown in Figure 52 is not too far from
reality. Since the dimensional reference length, the width of the lake, is
small compared to the radius of the earth, the centrifugal acceleration can be
neglected.
60
o
•H
P
u
o
ft
u
p
0)
o
p
Ch
o
i:3
O
•H
P
O
<D
W
I
m
o
o
p
m
a)
is
+3
,
ra
n
01
01
W
w
•r)
^
•
t)
H
•H
1
^
ir\
a)
(1)
,y
^H
01
g,
hi
•H
Ch
tq
o
61
T
D
1
EAST
Figure 52. Crosssection of a long, symmetrical trapezoidal
Lake Michigan model.
During the spring warming season, though the temperature structure in Lake
Michigan is in a continual changing semi steady state, the rate of heating in
creases with only some slight change in magnitude. Therefore the thermal body
force is a seasonal function varying slowly with time. In modelling the lake
with respect to the warmup season, all variables can be considered as time
independent without losing the good approximation.
The frictional forces are important near all boundaries. In the flow field
of a turbulent nature as in Lake Michigan, the frictional forces are expressed
as the space rate of change of the products of eddy viscosity and Reynolds
stresses in analogy to the laminar flow. The average horizontal eddy viscosity
in Lake Michigan is the order of 10^10^ cm^/sec and the vertical eddy vis
cosity (110^ cm^/sec) is much smaller. Eddy diffusivity may be of the same
order as eddy viscosity but it is a bit smaller in general. In our preliminary
model 10^ cm^/sec is used as the average value for eddy viscosity and eddy
diffusivity in both horizontal and vertical direction for ease in management.
In the future refinement, separate values for horizontal and vertical viscosity
are necessary in the Lake Michigan model. Eddy viscosity and eddy diffusivity
of the fluid depend not only on temperature (such as the molecular viscosity),
but change from place to place and time to time, depending on the local charac
ter and dimension of the flow. The overall average values of the order of mag
nitude of eddy viscosity and eddy diffusivity can be treated as "constant" for
application to the geophysical fluid problem (see section III of this report).
In general, though a weak thermocline exists during the later portion of
the warmingup season, the lake is homogeneous in this particular period. In
our preliminary model, we shall treat the fluid as homogeneous.
Boussinesq approximations are used in the formulation, and the equation of
state is a linear relation between the density and temperature as
p'  P^[l  a(T'  T^)],
62
where a. and T are specified in Table 51. Furthermore, the dissipation func
tion is neglected in the heat equation.
TABLE 5l» Nomenclatures and nondimensional parameters.
Symbol
Definition
X
Y
Z
q.
p
V
AAA
Ot
A
Aq
a
L
D
AT
Horizontal eastwest coordinate, meridional direction
Symmetrical north south coordinate, zonal direction
Vertical coordinate, negative is in the direction of gravity
Velocity with components (u^v,w)
Deviation of the pressure from hydrostatic pressure
Temperature
Constant lowest temperature near the bottom at the center of the
lake
Twodimensional gradient operator, i r + k r
ox oz
(N.B. The above symbols primed have dimensions, unprimed are dimen
Twodimensional Laplace operator,
(N.B. The
sionless. )
Unit vector in the (X,Y,Z) directions
Acceleration due to gravity
Coefficient of thermal expansion; always positive except below ^°C
Kinematic eddy viscosity
Thermometric eddy diffusivity
Angular velocity of the vertical component of earth ^s rotation in
Lake Michigan
Width of the lake^ reference horizontal length dimension
Depth of the lake, vertical dimension
Total horizontal temperature difference between the edge and the
central portion
65
TABLE 51 (Continued)
Symbol Definition
^ Stream function for the velocities
P Fluid density at temperature T
7 Depth to width ratio, 7 = D/L also nondimensional surface depth
\ Tangent of the slope angle with the vertical X = tan (see Figure
55)
C Dimension of velocity C = (XgDAT/2QL
p Nondimensional thermal Rossby number, p = agDAT/(2n)^ L^
G Nondimensional Ekman number, e = A/2nL^
Eddy Prandtl number, a = A/Aq
MATHEMATICAL FORMULATION
With all the considerations in the previous paragraphs, let us assume
that the fluid is constrained in a long symmetrical, trapezoidal, cross section
model of the lake. It is also subjected to two oppositely imposed horizontal
temperature gradients toward the center from both sides of the lake. The lake
center is taken as the center of the Cartesian coordinates (Figure 52). The
lake is rotating with angular velocity Q, (component of the earth's rotation nor
mal to the lake surface) about its vertical axis, which is antiparallel to the
gravitational force. The horizontal surface of the lake is considered free
with no wind stresses. The mean flow is assumed to be zonal, thus the basic
geostrophic thermal gradient relation is satisfied. All specifications of nomen
clatures are summarized in Table. 5l»
The governing equations of the lake model are then as follows:
q^.V^q* + 2nk X q* + i V^p»  ag(T^  T )k = AV»^J% (l)
AqV^^T'  q*V'T^ = 0, (2)
V»q^ = 0, (5)
6h
which are respectively the NavierStokes equations of motion, the energy equa
tion, and the continuity equation of the flow field. Note that the V and V*^
are twodimensional gradient and Laplace operators. The above equations are,
by axisymmetry, independent of Y. The boundary conditions are
T' = T , q' = at z' = 0, (k)
T' = (1  cos 2jtX')AT + T , ~^ = o at z = 7 (5)
2 o oz
i = , q'  at x= +1 L, (6)
which imply that a horizontal temperature profile of the form of an inverse
cosine curve is impressed upon the surface of the fluid of the lake, and the
lake is insulated at both edges with only a slight conduction between the lower
portion of the lake sides to the deeper layer of the earth. These differential
equations together with all boundary conditions are coupled high order and
nonlinear. Analytical direct mathematical solutions are very difficult. There
fore, some approximate approaches capable of making the equations more tractible
are preferred.
Nondimensional Equations
It is convenient in the approximate method to introduce nondimensional
variables and parameters. The following unprimed nondimensional variables are
defined:
rri t _ rp
L ^ D ^ AT ^ ^ C C^p ^
o
where C has the dimension of velocity. Since in our basic assumption the mean
flow field is in geostrophic thermal gradient equilibrium, then the geostrophic
thermal gradient relation:
20 — = ag —
dz "=" dx
leads to
agD
c =
^^^'Sx _ ggPAT
20, 2nL '
where D is the depth of the fluid.
65
Substitute all nondimensional variables into equations (l), (2), and (5)
with the change of vertical dimension D to the same dimension as the horizontal^
L, but using 2=7, the depth to length ratio^ to denote the free surface, the
governing equations become:
p(qV)q + k X q  kT + pVp = eS^q, (T)
gV^T  crpqVT = 0, (8)
Vq = 0, (9)
where
p _ ^ ^ is the nondimensional thermal Rossby number as the (h)
^ ^ of Hide^s (ISSk) parameter
e = p is the nondimensional Ekman number
a = A/Aq is the eddy Prandtl number
Both p and e are very small number and a has the order of unity. The fluid is
now bounded in
<x<7: ^0<z<7 where 7=7 ,
and boundary conditions are
q  , T = at z = 0, (lO)
■^0 ,T=(1 cos 2tix) at z = 7 (H)
oz 2
rVr 1
q = on the side boundary and ■T: = Oatx = +— . (12)
Since p is small, the ordinary perturbation method can be used. We may intro
duce a stream function by letting
66
t" =  JX Vt(xZ) + €"M y (13)
where the horizontal zonal velocity v is much greater than the meridional
velocity, u, and the vertical velocity, w, which will appear in the following
sections.
Expand all the variables as powers of p as
vn n (n)
^ n=0 "^ ^
n=0 ^
n=0 ^
where the superscript n denotes the nth order of p expansion. Equate terms in
the powers of p. The zeroth p order equations are:
W°^ +v^°^ T^°^ =0 (lU)
Z X
eV;°)  A°) =0, (15)
W°^ = 0. (16)
The first p order equations are
^4(1) (1) (1) r (O) (O) (O) (O) (O) (O) (O) (O),
Z X Z XZZ Z XXX X zzz X xxz
(17)
eW^) ^^'^  [/°)v(°^ ^(°V(°)] = 0, (18)
Z Z X Z Z
eW")  CT[^^°^T^°)  ^^^°^t(°)] = 0, (IV)
Z X X z
and similarly for the other high p order terms.
67
Boundary conditions are that all p orders of velocities, and in turn, all
porders of stream functions, should separately satisfy equations (lO), (ll),
and (12); and that all porders of temperature should satisfy (lO) and (12)
with (ll) for the zeroth p order and
T^^ = 0, (20)
for other higher ith porder temperatures, i = 1,2,3, etc.
SOLUTIONS OF THE ZEROTH p ORDER APPROXIMATION
Temperature
The zeroth porder temperature can be obtained from equation (l6), which
is simply the Laplace equation. The solution of equation (l6) satisfying boundary
conditions (lO), (ll), and (12) is
(o) Irz sinh 2nz ^  /^ v
r ^ = L  . ^ ^ cos 2TtxJ . (21)
2 7 smh 2tC7
The zeroth porder temperature is purely conduction.
Boundary Layer Scaling
Equations (ih) and (15) are coupled and show boundarylayer characteristics
because e is a very small number. In boundary layer analysis the boundary layer
thickness should first be determined as a function of the small parameter, e.
Equations (1^) and (15) imply that
eW°^ + tl°^ = . (22)
On the free surface and the rigid bottom we may scale the normal coordinates
near the boundaries as
cp :. €"^(z  9) (25)
where = on the bottom and = 7 on the surface.
68
Equation (25) substitutes into equation (22) with no stretching of the
horizontal coordinate. The small parameter €, multiplying the highly differ
entiated term is important only near the boundary where the frictional force is
the same order of magnitude as the other terms. This implies a = l/2.
On the inclined side boundaries ( see Carrier 1955)
where X is the slope of the slanted boundaries, and X = tan as shown in
Figure 52; let
I = ^""^(x + 0), (25)
and
where = +[2  X(7  z)] on each of the two side boundaries. By substituting
equation (25) into equation (22), to the lowest order of (e), we obtain b = l/2.
Therefore, all the boundaries have the boundary layer thickness of the order
Following the procedure of boundary layer additive type analysis, let
(o) (o) (o) (o) (o) (o) ...
r ^ = ^r ^ + ir ^ + ^r ^ + sr ^ + ^r . (26)
where the left subscript denotes the region. The jt is dominating in the
interior region and t^^^, i = 1,2,5^^ ^^^ significant only on their respective
boundaries. The various boundary layer regions and their widths are shown in
Figure 55' Coordinates for the regions of the lake cross section are summarized
in Table 52.
Interior Solutions
In the interior region, boundary effects are small. Equations (ik) and
(15), with the first terms neglected, become:
v(o) Jo) ^ ,(0) ^
^ ^ ^ = T^ ^ and ^\;^ ^ =
I X I z
69
JJL
(I)
12Y
(I) indicates the Interior region
(1) Indicates the surface region
(2) Indicates the bottom region
(3) Indicates the "eastern" boundary region
(4) Indicates the "western" boundary region
Figure 55. Boundary regimes of the Lake Michigan model,
TABLE 52. Local coordinates in different regimes.
Region
Coordinat
es
I
X
^
1
3£.
CPl = € (z •
■ 7)
2
X
1/2
cps = e z
5
ll
= e
.1/2
[X

1
2
+
Hy 
z)]
\=^
k
I2
= e
.1/2
[X
+
1
2

Hr 
z)]
Tl^Z
which yield:
and
v(o) l.cosh 27tz . >. n «/ \
E 2 smh 27t7 ^ ^ ^ ^
V°^ =B(x),
(27)
(28)
where A(x), B(x) are integration constants and possibly functions of x. A(x)
and B(x) can only be determined by matching up the interior solutions to the
solutions obtained from the boundary layers.
70
Boundary Layer Solutions
A. Region 1 — the upper boundary (free surface)
In the upper boundary region, the stream function and the zonal velocity
are,
f  jf + i¥ ^ V  V + iv
and the vertical coordinate is stretched as shown in Table 52. Equation (22),
using the local coordinates, yields
Since it^ ^ vanishes away from the upper boundary, no nonzero forcing term can
be satisfied after integration of equation (29). Therefore equation (29) can be
simplified as:
Together with the appropriate boundary conditions, as
5.^°^ ,(o) ^i^°^ , ,
or equivalently,
equation (30) can be solved and we obtain:
(o) cpi/f2 . , I , ,
i\i/^ '^^ = GTce^^'^ sin 2Ttx cos cpi/v2 , (32)
and
( o) 1/2 Tte . ^ / cpi cpix
iv^ = eV^ ^^ sin 2^.x(cos ^ + sin ^)
71
(35)
B. Region 2— the lower boundary (rigid bottom)
Stretching the vertical coordinate^ equation (22) leads to an equation
similar to (50). With the lower boundary conditions
(o) (o) S^l;^''^
"^ ^""t = at cps = , (5U)
we obtain
2t' =€:te ^^'^ sm 2Trx(cos ^ + sin ^) , (55)
and
M°^ = cV^vfl «e^^/^ sin 2«x cos ^ .
(56)
Matching the boundary solutions with that of the interior, the boundary
conditions in equations (5I) and (^k) also yield the two integration constants
in equations (27) and (28) as
A(x) = sin 2Ttx(€V2 fz   — ,^ ^ ) ,
^ ^ ' 2 sinh 2Tr7^ ^
and
B(x) = €tt sin 2Ttx .
Then the interior solutions are:
(o)
j^ = Gjr sin 2:n:x , (57)
(o) , 1/2 nr cosh 2jtz  1 , . ^ ...
_v' = (g / \r2 +• — . , ^ ) sin 2nx . (58)
I ^2 smh 21(7
C. Regions 5 and h — eastern and western rigid side boundaries
Using the local stretched coordinates normal to the side boundary, equation
(22), to the lower order of €, leads to
72
Oil
in region 5»
With the boundary conditions,
So)
( o) St , , .
v^ ^ = ^— = at li = (UO)
051
we obtain.
,(°) e(l + A.^)k kii/y/i" . ^ ,1 ,, ... cosh 2Ttnj  1 r^n/ Mi • Hi\
(o) ki/>/2' ki . ^ a , . ,,cosh2itTii 1 r?r 1/2,
where k = \ /(I + X ) . By symmetry, similar solutions are obtained in region
k as
4'
(o) e(l + \^)k )^iz/{2 . ^ a ,, , , , cosh2rtT]2l r~,. kfo . kp.
,(o)_ k2/V2"_„ k
4V • = e ' ■ cos
(111.)
All the boundary layer solutions are significant only on their respective
local boundaries and go to zero away from their own boundaries. From equation
(26), the general solution of the zeroth p order approximation is:
75
.^^°^ = e[n sin 2^[l . e^'^/^^ cos ^^ . e"^/'^ (cos ;#= . sin ^) ]
_.[..l/2X(..)„ ^ (^^,
J°) • o • r cosh 2itz  1 ^ 1/2 ,r rt J^e^ '^^ rrr z/\f2e z ,
^ = ^^^ "'^f 2 sinh 2.r """^ ^^W^^ ^^^ °°^;??]
^ {sin 2n[i  A.(7  ,)](<^°^^2^^ 1 + ^1/2^), ^^k[x4l/2X(rz)]/v[2r
2 2 sinh 2tC7
 cos kfx + 1/2  X(y  z) 1 k[xl/2+\(rz) ]/\[2r kfx  l/2 + A.(y  z) 1
cos ^  e cos ^ ) '
and equation (21) is the zeroth porder temperature.
RESULTS OF THE ZEROTH p ORDER SOLUTION
For a small thermal Rossby number^ p^ and a high Taylor number^ equivalent
to € ^ our solutions are within the range of lower symmetrical regime in a
horizontal differentially heated rotating fluid (Fowlis and Hide I965). The
governing equations are linearized by ordinary, perturbation expansion with p.
Solutions are obtained by further singular perturbation expansion for each
order of p in performing boundary layer analysis with respect to the Ekman
number €.
Temperature
The zeroth porder temperature in equation (21) is obtained from heat con
duction, neglecting the convectional and the frictional effects. The temperature
distribution is a linear function of depth peaked on the surface edge with some
horizontal variation minimized at the central portion of the lake. Figure 5^4
presents the horizontal temperature distributions of available data at the Great
Lakes Research Division from 1965 to 1968 during the warmup season. Figure
55 shows the vertical temperature distributions during the same period of time
in the central portion of Lake Michigan. Generally speaking, both the horizontal
and the vertical temperature distributions show the same pattern as the solution
indicated. The solution is in good agreement with the mean synoptic temperature
measurement, especially in the interior region. On the boundaries, either hori
zontal or vertical, the agreement is not as good as the interior. On the upper
::::!::9 12.0 10.0 10.0
43** 30'
JUNE 1965
GRAND HAVEN
43«00*
19.0
17.0
 43*^30'
JULY
GRAND HAVEN
43° 00
■ 43»30'
AUGUST
\v
\ 1/'^^^^^^^ HAVEN
' ' \ V yy^y^y 43*»oo'
(a)
Figure '^k. Surface temperature patterns during spring warming
of Lake Michigan from available data, I95U68.
75
87O00'
4.0
43*»30'
MAY 1966
GRAND HAVEN
::•::.' 43« oo'
8.0 .6.0
V.
^^^
_ — —
/
7.0
\
/
/
/
/
1
\
7.0
\
/
\
43** 30'
JUNE
GRAND HAVEN
. 43«00'
43'* 30'
AUGUST
GRAND HAVEN
X 43** 00'
20.0
(b)
Figure 5^ (Continued),
16
■ 43° 30'
MAY 1967
GRAND HAVEN
43° 00'
■ 43° 30'
mi I '
f / \ I \ \wmm JUNE
il ' i 1 1 W
X'>7 / . / /I5 0)^„
I }i50^^^^^^^ HAVEN
Figure 51^ (Continued).
77
6.0 4.0
► •"•"•'•"• «
•••••«
• ••••••• aV ■
x:::v:vy * 14.0 \
• •••••• • \ ' *\j,\j %
12.0
87*»00'
iig^im::.. MAY 1968
43*' 30'
GRAND HAVEN
43** 00'
43** 30'
Figure ^k (Concluded).
78
CO
o
CO
LiJ
UJ
<
a
o
•H
P
^
O
Ph
H
(d
^
•P
S
<u
C)
•
00
«H VO
o
1
\r\
0)
M3
^
ON
H
4:>
o
•\
g
t
p
•H
QQ
Jh
(U
ft
b
bO
P
d
cd
^
g
0)
i
Q)
+3
bO
a
H
•H
cd
O
1
•H
TJ
1^
0)
§
>
ttO
:d
k
o
g
sj a; 9UJ
79
CD
Z
o
CD
LiJ
CO
Z)
llJ
<
80
MILWAUKEE
MUSKEGON
MAY
MUSKEGON
16° I7*»_
IP 17°
JUNE
19° 18*
1967
I9« 19°
2»l8°
AUGUST
Figure 55 (Continued).
81
00
CD
CD
o
UJ
o
00
o
g
5
^
J^
UJ
1
^^
1 1
S
SJ949UJ
82
horizontal boundary, the real temperature decrease with depth is more rapid than
the relation predicted by the solution. On the lower bottom layer, the real
temperature decrease is much less. On the two slant side boundaries, real
temperature gradients exist quite close to the shore, and in the central portion
of the lake the temperature change is very small. The temperature variation is
much greater than the solution yields. All these observations may indicate:
(1) the impressed surface temperature is not an exact inverse cosine profile,
and (2) the temperature is modified by convection as a result of the coupling
of the temperature with the inertial terms in the equation of motion, particularly
on the boundary layers. For mathematical simplicity, by smoothing the tempera
ture gradient, the inverse cosine temperature distribution is not too far from
reality. As for the importance of convection in determining the temperature,
higher porder approximations should be attempted. Since the zeroth .porder
velocities have boundary layer characteristics, and these velocities together
with the zeroth porder temperature gradients are forcing the first porder
temperature, it is expected that the higher porder temperature will show
boundarylayer characteristics and modification of the overall temperature dis
tributions.
Velocity
From the basic realistic assumption of the geostrophic thermal gradient
balance, neglecting variable surface wind stresses, the zonal velocity dominates
the overall circulation of the whole lake as expected. In the interior region,
the pressure force balances the Coriolis force due to the zonal velocity together
with the thermally induced body force. This basic assumption of geostrophic
thermal gradient relationship (the balance of the vertical zonal velocity gradient
and the horizontal temperature gradient) must be maintained. The interior solu
tion of the zonal velocity is shown in equation (38). It has the maximum magni
tude around onefourth the width of the lake from the shore edge near the sur
face and approaches minimum as x approaches zero and half indicating the middle
portion and both edges of the lake. During this warmup season, on the eastern
portion of the lake where x is positive, there is a northward zonal current.
On the western portion of the lake, where x is negative in the horizontal coor
dinate of the lake model, there is a southward current. Both the horizontal and
the vertical velocities are much smaller (in the order of e or ei/s) than the
zonal velocity. From the interior stream function equation (57), there is no
horizontal velocity (u = SVSz) in the interior region. The vertical velocity,
w =  d\f/Sx, leads to a broad but very small vertical circulation in magnitude
in the interior region given as
(o) P
^w' ' = € 2Tr cos 2Trx ( J4.7)
with a relatively strong sinking in the middle of the lake, and slight updwell
ing along the edges.
85
In region 1, that is^ in the free surface layer^ the zonal velocity is
strengthened by adding to the interior zonal velocity a boundary zonal velocity
of order (g*/^) in magnitude. The direction of the general zonal flow is the
same as that of the interior, which shows a cyclonic circulation. The
meridional velocity in the upper boundary is different from the interior where
no eastwest horizontal component develops. On the surface layer^ there is a
meridional velocity
..(o) St 1/2 ^ 1/2 ^ ^i/{^ . ^ / 9i . ^i\
^^' = ^ = € /^ L_ = e/^:^ e^^/^ sin 2^(cos ^  sm ^)
(^8)
which flows from shore boundaries due east and due west^ toward the central
portion of the lake. The vertical velocity in region 1 is much smaller than
the horizontal velocity as the continuity equation implies. The vertical velo
city is:
(o) (o) (o) (o) ^2 ^1/^2 cpi ^ .,^.
W^ ^ = W^* + iW^ I^ ^ + G 2Tr e^^^"* cos ^ cos 27tK (1+9)
which is the continuation of the interior velocity with only a slight increase
in intensity.
In region 2, the bottom layer, there is a meridional flow of an equal or
der of magnitude but opposite in direction to that of the surface, with a ver
tical velocity similar to that of the surface. The zonal velocity in the bottom
layer is much smaller, as expected. It can be noticed that the residual zonal
velocity near the bottom is the difference of the interior and the boundary
zonal velocities.
In the side boundaries, regions 5 sind U, the local coordinates (^r) are
stretched normal to the inclined surfaces. From equation (^l) the meridional
velocity in region 5, to the order of e, is:
/ 1/2 . S S V ,( o)
,(0) _ ^t^°^
that is
(0)
eV^ k^(i + x2)^ki.
. k^i/>/2' . k^i . ^ , v.cosh 21^111 ^ rT\
)e sxn f sxn 2^7 " ri3)( , ^^^^\^  <^)
(50)
&k
This meridional current flows away from the shore with a moderately high speed
in the side boundary as i is negative when X is positive, and vanishes in the
interior region. By symmetry, a similar flow away from the west coast toward
the center of the lake exists in region k.
The local vertical velocity in the eastern alongshore region, to the order
of e, is
,w
(o)
ax
= ._ei/2 k^d + x^)
e^^^/ ^ sin^ sin 2TtX(r
Tll)(
cosh 2Trr]T  1
2 sinh 2ny
(51)
which represents sinking in the side bounda ry layer . The local vertical velo
city in the alongshore boundary region is opposite to the vertical velocity in
the interior. Therefore in the alongshore boundary region^ the overall vertical
velocity may be small. Similar symmetrical phenomena exist in region k.
The local zonal velocity in regions 5 and k, as shown in equations (^2)
and {^h), are opposite in direction to the general pattern of the zonal circu
lation extended to the boundary from the interior. In the east and west along
shore boundary region^ the current direction is dependent on the difference in
magnitude of the general zonal current and the local zonal current. Therefore^
the overall zonal currents in these two alongshore boundary regions may flow
either southward or northward. Theoretically^ boundary currents very close to
the shore should be smaller.
In summing up the foregoing current analyses^ the general horizontal zonal
circulation is shown in Figure 56. The meridional and the vertical circulations
are shown in Figure 57a When the overall temperature in the lake is less than
i^.°C^. the coefficient of thermal expansion of water^ a, is negative. Hence the
meridional and vertical circulation patterns of the lake are changed as approx
imately shown in Figure 5Tt>*
WEST
SHORE
A^\
I
r i
J I I
(1)
(I) ♦
(2)
f f.l
1(3)
EAST
SHORE
Figure 56. General horizontal zonal currents.
85
WEST
:^ * (i)* * 5= 5rT=2
i i i(I)l I i f
I t
*  ^ *. *(2) * _^ _^ ^ ^
v. V> V^ v^ v^ v^ v^ *^ ^^ V^^ V
a) . For a positive
EAST
^ ^ * — —. (1) * 
WEST
i t »
t I'^'Im
^ ^ ^ ^ (2)
. . V^ V. v^ V.
b) • For a negative
Figure 57. Meridional and vertical circulation.
EAST
Since we have not analyzed in detail the magnitude of each component of
velocity. Figures 56 and 57 only roughly show the possible general circula
tion pattern. A synoptic surface current computed from dynamic heights in
June 1955, reproduced from Ayers et al. (1958), is shown in Figure 58. At
the central portion of the lake, the general circulation pattern in Figure 56
is in good agreement with that of Figure 58. A series of drogue studies for
surface current measurement are shown in Figure 59 for the eastern shore area,
and in Figure 510 for the western shore area. Tables 55 ^^^ 5^ give the
elapsed time for the drogue runs and the average drogue speeds for the east
shore and westshore experiments, respectively. Generally, drogues drift toward
the north near the eastern shore and toward the south near the western shore,
which is in good agreement with the theoretical solutions. Winds are also in
dicated in Figures 59 and 510. In some cases, currents were flowing directly
against the winds. In studying the thermal bar phenomena (the temperature gra
dient zone) in the Great Lakes, Rodgers (1965) indicated that the circulation
associated with the thermal bar should be in the pattern as shown in Figure 511
(Figure 7 of Rodgers* paper). When compared to our theoretical current pattern
in a vertical section shown in Figure 57ti the major characters of these features
agree well in general. Of course, our model does imply that this kind of circu
lation pattern exists only when the temperature in the thermal bar is lower than
4°C. It is expected that when the sharp temperature gradient exists around the
k°C zone, this kind of circulation pattern is more augmented by the acceleration
of the vertical velocity in the boundary by the negatively buoyant force.
86
Figure 58. Surface currents in Lake Michigan, 29 June 1955
(after Ayers et al. 1958)
87
86°30'
43°10'
(a)
/
\
B.
86^20
(b)
86°30'
43**00'
/
K. L.
M.
V
86°20»
43^00*
Figure 59 • Drogue studies for surface current mea,surements near
eastern shore, (a) 2 May I968 (wind N, 10 kts?). (b) 2k June
1968 (wind SW 15 kts).
88
(c)
66°30»
43°10»
Date
TABLE 53 Drogue Studies on the East Side.
Wind Drogue Station Start Finish cm/sec.
2 May '68 N10
24 June
SW15
A.
B.
C.
D.
E.
F.
G.
H.
I.
J.
K.
L.
0950
1003
1009
1015
1020
1026
1031
1032
1027
1022
1017
1012
1504
1514
1522
1534
1545
1554
1559
1502
1520
1548
1602
1612
7.14
5.80
5.80
3.12
6.24
6.24
4.01
8.94
5.80
7.70
10.26
8.96
30 July SE20
M.
1007
1623
9.41
N.
1002
1637
9.41
0.
1005
1328
28.9
P.
0952
1300
19.5
Q.
Lost
R.
0931
1220
10.1
S.
0918
1235
12.7
(c) 50 July 1968 (wind SxE, 20 kts?)
Figure 59 (Concluded).
89
43*'00
Milwouke*
(a)
CALM
E.
F.
T^
^
TABLE 54 Drogue Stations on the West Side.
Wind Drogue station Start Finish cm/sec .
10 May '68 Calm
Miles
27 June
87^50'
NNW6
43*^00'
\
A.
0800
1506
3.12
B.
0815
1520
4.73
C.
0821
15 32
8.02
D.
0837
1548
4.91
E.
0848
1604
7.58
F.
0859
1620
5.36
G.
0912
1635
4.01
H.
1245
1848
4.91
I.
1251
1820
9.86
J.
1301
1813
10.74
K.
1308
1805
7.14
L.
1314
1856
4.91
M.
1321
1840
8.47
N.
1327
1831
7.14
(b)
43°
Figure 510. Drogue studies for siirface current measurements near
western shore, (a) 10 May I968 (wind slight), (b) 27 June I968
(wind NxW, 6 kts).
90
WARM^N^
Figure 511. Circulation associated with "Thermal Bar"
(after Rodgers I965).
THE FIRST p ORDER SOLUTIONS
As we have already mentioned, the first p order governing equations of the
model are equations (1?)^ (18)^ and (19). The forcing terms which are the
products of the zeroth porder temperature gradients and the velocities (as
shown in equation (19))^ are of the order unity in the interior region, of
the order (e"""/^) for the individual local boundary layer regions at their re
spective boundaries, and with another forcing term of the order unity due to
the overall interior influence. Since the forcing terms are of the boundary
layer type, it is expected that the first p order temperature will have strong
boundary layer characteristics. In each of the boundary regions, the first p
order temperature can be separated into two parts; one part is defined by the
overall existing zeroth order interior stream function (^t^^^), and the other
part is due to the local stream function which is important only in its own
boundary region. In this way the boundary layer portion which is forced by
the terms of order of (e ^/^) can be treated with boundary layer analysis.
After applying the respective local boundary conditions, solutions are obtained
for all the first p order temperatures influenced by local zeroth p order stream
functions. Then at all the new boundary conditions, solutions are obtained by
applying the first porder solution which has the boundary layer characteristic
as a perturbation of that part of the temperature under the Laplace operator
(which has no boundary layer characteristic). In the interior region the tem
perature is separated as a general solution and a particular solution, with the
forcing terms of the order of unity operating only on the particular solution.
After transferring the nonhomogeneous terms as a particular solution into the
new boundary condition and matching with the boundary solutions to the order of
unity, the overall first porder temperature can be solved. The result of the
first porder temperature turns out to be of the order of unity and there are
some modified variations of the order ( e^/^) in the boundary regions.
91
After obtaining the first porder temperature, the first porder stream
function and the zonal velocity can be obtained from the coupled equations (I7)
and (18) following the method used in the zeroth porder solutions, and with
more complexity.
CONCLUSIONS
After the winter mixing, water in Lake Michigan is almost homogeneous.
The spring warming season begins during late April or early May, ' Beginning at
this time the lake starts to warm from the shallow water along the shore
progressively outward to the main body of water which is still at k^'C or lower.
During this period, the land and air are much warmer than the water. A stable
situation usually exists in the atmosphere above the lake. It is expected that
much less energy is transferred from the wind to the lake water. A sharp tem
perature gradient always exists within a few miles along the shore which leads
to a balanced geostrophic thermal gradient model of Lake Michigan in a simpli
fied long trapezoidal shape. A double perturbation method is used in obtaining
approximate solutions from these coupled nonlinear momentum and energy equations
by expanding first with the small thermal Rossby number, p, and then, to each
order of p using singular perturbation, with respect to the small Ekman number,
€. All boundaries are found to have the same order of boundary layer thickness
in the order of ( € ^f). The zeroth porder temperature is found to be due en
tirely to conduction. The first porder temperature shows that convection is
important and it has boundary layer variations. All the velocities do have
boundary layer characteristics. The general zeroth porder horizontal zonal
velocity shows a cyclonic circulation pattern with variations in the along
shore boundaries. When the overall temperature of the lake is still lower than
k°C, the horizontal meridional velocities and the vertical velocities show a
broad upwelling of a small order of magnitude in the middle portion of the lake
and a sinking zone along the edge of alongshore boundary layers coupled with the
divergent horizontal surface meridional current and the convergent bottom cur
rent. This forms two gigantic but weak counterrotating cellular motions. Two
horizontal meridional currents flowing outward from their respective shores in
the two alongshore boundary regions to joint the thermal bar convergences sug
gest the possible existence of another upwelling very close to the shore bot
toms. Though detailed characteristics have not been wholly demonstrated, the
general circulation pattern compared with results obtained by Ayers et al.
(1958) sind Rodgers (1965) indicate very good agreement. From a physical point
of view, in the spring warming season the shallow water alongshore region re
ceives a great amount of heat energy from the air and the land. The temperature
in these regions increases very rapidly while that of the main body of water
does not respond as quickly and remains at the temperature of maximum density.
This temperature gradient induces an inclined surface water slope tilted up
toward the land, which induces a "Hadley cell" (1755) The thermal current is
originally flowing outward from shore but is deflected. by the Coriolis force
which results the northward or southward zonal current. Due to the relatively
small dimensions, the boundary effects are more important in lakes than in
92
oceans. Though we have not analyzed the numerical magnitude of all variables^
the general circulation pattern has demonstrated all the main expected results,
FUTURE PIAN
Our preliminary model is too rough to obtain detailed characteristics of
the lake circulation either in magnitude or in flow pattern. It is my hope to
be able to further refine the model both in dimension and magnitude of some var
iables in similation of the prototype of the lake. For example^ the horizontal
dimension and the vertical dimension should be scaled separately; the vertical
velocity should have a different velocity scale than the horizontal one; the
horizontal eddy viscosity and the vertical eddy viscosity should have different
mean values. As in perturbation procedures, higher orders of approximation are
very important. The detailed analyses of the first gorder solutions and the
estimation of approximate error from the second Porder equations are scheduled
to be carried out. The combined final solutions should be programmed into com
puter calculations and the results should be compared with all available field
data including some of the current meter buoy stations of the Great Lakes
Illinois River Basins project in 196^. Furthermore, I also intend to consider
the thermal body force as a slowly time dependent function and to extend the
application of the model for the entire year.
95
References
AYERS, J* C D. E. CHAHDLER, G. H. LAUFP^ C. F. POWERS, and E. B. HENSEN.
1958. Currents and water masses of Lake Michigan. University of Michigan,
Great Lakes Res. Inst. Pub. No« 5, I69 p.
BELLAIRE^ F. R. I965. The modification of warm air moving over cold water.
Proc. 8th Conf., Great Lakes Res., University of Michigan, Great Lakes Res.
Div. Pub. No. 15, pp. 2^19256.
CARRIER, G. F. 1955* Boundary layer problems in applied mechanics. Adv. in
Appl. Mech., 3:119.
FOWLIS, W. W. and R. HIDE. 1965. Thermal convection in a rotating annulus of
liquid: effect of viscosity on the transition between axi symmetric and non
axisymmetric flow regime. J. Atmos. Sci. 22:5Ul.
HADLEY, G. 1755. Concerning the cause of the general tradewinds. Philo
Trans. Roy. Soc.j, London, 59:5865.
HIDE, R. 196U. Experiments on baroclinic instability and "vacillation": a
review. Trans. Am^ Geophys. Union ^5: 1^+21^+5*
NOBLE, V. E. and R. ANDERSON. I968. Temperature and current structure in the
Grand Haven, Michigan, vicinity during thermal bar conditions. In press, Proc.
11th Conf. on Great Lakes Res.^, Internat. Assoc, for Great Lakes Res.
NOBLE, V. E. 1967. Evidences of geostrophically defined circulation in Lake
Michigan. Proc. 10th Conf. Great Lakes Res., Internat. Assoc, for Great Lakes
Res., po 289298.
RODGERS, G. R. 1965. The thermal bar on the laurentian Great Lakes. Proc.
8th Conf. on Great Lakes Res., University of Michigan Great Lakes Res, Div.
Pub. No. 15:558565.
_____^ . 1966. The thermal bar in Lake Ontario, spring 1965 and winter
196566. Proc. 9th Conf. on Great Lakes Res., University of Michigan, Great
Lakes Res. Div. Pub. No. 15:56957^
STOERMER, E. F. I968. Near shore phytoplankton populations in the Grand Haven,
Michigan, vicinity during thermal bar conditions. Proc. 11th Conf. Great Lakes
Res., Internat. Assoc. Great Lakes Res. In press.
9h